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General  note  1: 


Vectors  and  dyadics,  which  are  independent  of  reference  frame,  are  denoted  by  a  bold  char¬ 
acter  and  an  arrow  (e.g.  h  and  I).  Unit  vectors  are  similarly  denoted  by  bold  characters 
and  hats  (e.g.  bi).  Matrices  representing  vectors  and  dyadics  in  a  particular  reference 
frame  are  denoted  by  the  same  bold  character  without  the  arrow  or  hat.  The  particular 
associated  reference  frame  will  be  stated  in  the  accompanying  text  or  implied. 


General  note  2: 

A  tilde  superscript  on  a  vector  or  a  matrix  representing  a  vector  (eg.  M)  indicates  that 
it  is  the  dimensional  form.  Accordingly,  vectors  and  matrices  representing  vectors  without 
tilde  superscripts  can  be  assumed  to  be  non-dimensionalized  unless  specifically  indicated. 
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Abstract 

Attitude  maneuvering  of  three  axis  stabilized  spacecraft  is  typically  accomplished 
using  a  combination  of  expendable  propellant  thrusters,  magnetic  coils,  momentum  storage 
mechanisms  such  as  momentum  wheels,  reaction  wheels,  or  control  moment  gyros,  or 
passive  methods  such  as  gravity  gradient  stabilization. 

With  a  few  exceptions  (most  notably  the  Hubble  Space  Telescope),  most  three  axis 
stabilized  spacecraft  in  service  today  still  rely  heavily  upon  expendable  propellant  thrusters 
for  primary  attitude  control,  using  momentum  storage  devices  only  to  offset  small  pertur¬ 
bative  torques. 

While  thrusters  allow  rapid  changes  in  attitude,  their  reliance  on  limited  propellant 
supplies  can  result  in  shortened  useful  lifetimes  for  spacecraft  employing  them  exclusively. 
Passive  attitude  control  methods  such  as  gravity  gradient  stabilization  or  utilization  of 
solar  pressure  on  solar  arrays  tend  to  be  slow  and  are  dependent  upon  interaction  with 
uncontrollable  forces.  Maneuvering  using  magnetic  coils  is  likewise  somewhat  slow  and 
becomes  less  effective  with  altitude. 

Momentum  transfer  devices,  on  the  other  hand,  rely  upon  a  reusable  energy  source  - 
electrical  power  -  provided  by  solar  arrays  and  rechargeable  batteries,  and  can  be  used  to 
effect  large  angle  reorientations  relatively  quickly.  However,  due  to  the  inherent  gyroscopic 
coupling  for  a  rigid  body  with  embedded  momentum  wheels  and  limited  control  authority 
in  terms  of  allowable  torques  on  the  wheels  and  maximum  allowable  wheel  speeds,  im¬ 
plementation  of  “optimal”  control  laws  for  these  devices  can  be  somewhat  complicated. 
Calculation  of  wheel  torques  for  large  angle  maneuvers  based  on  eigen-axis  rotations  or 
other  currently  used  control  laws  are  beyond  the  on-board  processing  capability  of  most 
current  satellites.  This  results  in  the  vehicles  being  dependent  on  the  ground  segment 
computers  and  ground  commanding  for  such  maneuvers. 

This  thesis  addresses  sub-optimal  employment  of  3  momentum  wheels  for  large  angle 
reorientation  of  rigid  spacecraft  with  minimal  induced  spacecraft  motion  during  maneuvers. 
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In  addition  to  development  of  general  theory  for  3  wheel  vehicles,  simulation  results  for 
a  vehicle  using  momentum  wheels  for  secondary  attitude  control  (GPS  Block  HR)  are 
compared  to  results  for  a  vehicle  using  them  for  primary  attitude  control  (the  Hubble 
Space  Telescope),  to  demonstrate  practical  applications  and  limitations. 

While  the  control  laws  were  developed  assuming  no  external  perturbing  torques  on 
the  vehicle,  reorientation  scenarios  were  run  both  in  a  torque  free  environment  as  well 
as  an  environment  with  simulated  gravity  gradient  and  solar  pressure  torques.  The  goal 
was  primarily  to  show  the  growth  of  vehicle  angular  velocities  and  again  demonstrate 
limitations  of  the  derived  control  laws. 

The  results  indicate  that  for  real  spacecraft  with  limited  wheel  momentum  storage 
capacities,  there  is  a  significant  trade-off  between  maneuver  times  and  required  wheel 
torques,  and  that  final  state  errors  (angular  velocities)  increase  with  increasing  wheel 
torques.  Nonetheless,  the  simulations  demonstrated  that  large  angle  maneuvers  can  be 
performed  for  both  GPS  Block  HR  and  Hubble  Space  Telescope  in  reasonable  times  and 
with  small  angular  velocities  using  the  sub-optimal  control  law.  However,  gravity  gradient 
and  solar  pressure  torques  tended  to  cause  larger  fluctuations  in  total  angular  momentum, 
angular  velocities,  and  final  state  errors  for  the  Hubble  Space  Telescope. 

HopefuUy  the  results  of  this  report  can  be  used  as  a  basis  for  future  analysis  of  this  or 
similar  sub-optimal  control  approaches.  Recommendations  for  follow-on  research  include 
comparison  of  this  approach  to  optimal  control  formulations,  as  well  as  consideration  of 
vehicle  kinematics  in  the  control  law  formulation  to  reduce  final  attitude  errors.  This  can, 
in  turn,  be  used  in  the  development  of  closed  loop  controllers  that  will  minimize  deviations 
from  the  sub-optimal  control  trajectories  in  the  presence  of  external  perturbations. 
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SUB-OPTIMAL  CONTROL 


OF 

RIGID  SPACECRAFT  REORIENTATION 
USING 

THREE  MOMENTUM  WHEELS 

1.  Introduction 

1.1  Background 

Spacecraft  attitude  dynamics,  stability,  and  control  has  obviously  been  an  area  of 
avid  research  since  the  launch  of  the  first  artificial  satellites  in  the  late  1950’s.  Although 
early  satellites  had  no  need  for  accurate  attitude  control,  the  rapid  advancement  in  space- 
based  communications  and  remote  sensing  applications  has  resulted  in  the  need  for  robust 
and  efficient  attitude  control  systems  in  modern  spacecraft. 

Most  communications  and  remote  sensing  spacecraft  in  operation  and  in  development 
today  employ  three-axis  stabilization  due  to  mission  mandated  pointing  accuracy  require¬ 
ments  and  the  associated  need  to  minimize  vibrations  in  large  flexible  appendages.  With 
a  few  notable  exceptions  (e.g.  Hubble  Space  Telescope),  three-axis  stabilization  schemes 
typically  utilize  momentum  storage  devices  such  as  reaction  wheels,  momentum  wheels, 
or  control  moment  gyros  only  to  offset  the  long  term  effects  of  small  external  perturbing 
forces,  and  accomplish  large  angle  reorientations  using  expendable  propellant  thrusters. 
While  certainly  efficient  in  terms  of  controllability  and  time  required  to  complete  maneu¬ 
vers,  the  use  of  thrusters  limits  the  effective  life  of  satellites,  results  in  high  maneuver  costs, 
can  induce  undesirable  motion  in  spacecraft  appendages,  and  can  in  some  cases  adversely 
affect  the  operation  of  sensitive  on-board  systems  (like  optical  telescopes  in  the  case  of 
Hubble).  It  is  thus  in  general  cost  effective  and  in  some  cases  mission  critical  to  use  mo¬ 
mentum  storage  devices  not  only  for  fine  attitude  control  but  also  for  primary  attitude 
maneuvering  whenever  practical. 
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Many  of  the  analytic  solutions  to  spacecraft  attitude  control  problems  to  date  have 
focused  on  the  relatively  simple  case  of  the  gyrostat,  which  is  typically  used  to  refer  to 
“dual-spin”  spacecraft  composed  of  an  inertially  fixed  segment  (platform)  and  a  spin¬ 
ning  segment  (rotor)  which  provides  gyroscopic  stiffness.  Several  textbooks  (most  notably 
Hughes  [11])  journal  papers  (including  Hall  [8], [9],  Hall  and  Rand  [10]),  and  AFIT  theses 
(Kinney  [14],  Kowall  [15],  and  Tsui  [22])  have  addressed  exact  or  approximate  analytic 
solutions  and  related  stability  for  both  single  and  two  rotor  gyrostats.  Many  researchers 
have  also  addressed  developing  control  laws  for  multiple  momentum  wheels  from  an  op¬ 
timal  control  perspective  (a  significant  contribution  is  provided  by  Junkins  and  Turner 
[12]).  Unfortunately,  the  resulting  calculations  can  be  computationally  intense,  and  thus 
relegated  to  ground  segment  computers  for  practical  application.  This  increases  satellite 
dependence  on  the  ground  segment,  limits  their  autonomy  (and  thus  survivability  and 
flexibility),  and  increases  operations  costs  for  these  vehicles. 

However,  Hall  [8]  did  provide  the  framework  for  development  of  a  class  of  sub-optimal 
maneuvers  based  on  his  analysis  of  single  and  two  rotor  gyrostats.  This  stationary  platform 
maneuver  is  a  rest-to-rest  maneuver  (inertially  fixed  at  beginning  and  end)  during  which 
angular  velocities  of  the  platform  are  kept  small  throughout  the  maneuver.  This  has 
the  advantage  of  limiting  induced  vibrations  in  flexible  appendages,  and  is  also  easily 
calculable.  Thus,  such  a  control  law  could  easily  be  incorporated  into  the  attitude  control 
logic  of  future  spacecraft  and  calculated  on-board,  thereby  increasing  autonomy  of  these 
vehicles  as  well  as  the  mission  life  by  using  momentum  wheels  for  large  angle  maneuvers 
instead  of  thrusters. 
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1.2  Goals  of  This  Thesis 

The  primary  objectives  of  this  report  are  three-fold:  To  extend  Hall’s  research  by 
developing  a  “simple”  sub-optimal  control  law  for  spacecraft  using  three  momentum  wheels, 
to  assess  the  utility  and  limitations  of  this  control  law  for  “real”  spacecraft  in  a  torque  free 
environment,  and  finally  to  qualitatively  address  the  additional  limitations  of  the  control 
law  imposed  by  the  addition  of  external  perturbing  torques  such  as  gravity  gradient  and 
solar  pressure. 

This  work  will  thus  hopefully  serve  as  a  foundation  for  future  research  using  this 
or  similar  control  approaches  by  identifying  the  strengths  and  weaknesses  as  currently 
implemented  and  providing  specific  recommendations  for  improvement. 

1.3  Methodology 

A  computer  simulation  program  was  developed  by  the  author  for  use  in  this  research 
effort.  The  code  is  included  in  Appendix  D.  This  program  calculates  the  physical  proper¬ 
ties  for  the  GPS  Block  HR  and  Hubble  Space  Telescope  vehicles,  sets  initial  orbital  and 
simulation  parameters,  numerically  integrates  the  equations  of  motion  for  the  vehicle  (per¬ 
turbed  or  unperturbed),  and  provides  both  textual  and  graphical  output  of  the  results. 
Data  from  the  simulation  runs  are  the  main  subject  of  analysis  in  this  report. 

Chapter  2  of  this  report  focuses  on  the  derivation  of  the  equations  of  motion  of  an 
Earth  orbiting  rigid  spacecraft  with  fixed  momentum  wheels,  in  the  presence  of  gravity 
gradient  and  solar  pressure  torques.  It  is  these  equations  which  are  integrated  by  the 
aforementioned  simulation  program.  While  mainly  a  compilation  of  results  from  previous 
research,  the  chapter  provides  a  coherent  and  logical  summary  of  the  pertinent  equations, 
thereby  reducing  the  need  for  the  reader  to  refer  to  outside  sources. 

In  Chapter  3,  Hall’s  development  of  sub-optimal  control  laws  for  a  spacecraft  using 
two  momentum  wheels  is  extended  to  a  spacecraft  using  3  (or  more)  wheels,  resulting  in 
a  wheel  torque  control  law  that  results  in  “smooth”  large  angle  reorientations  from  one 
inertially  fixed  state  to  another  with  only  small  induced  body  angular  velocities. 
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Chapter  4  compares  this  sub-optimal  control  law  to  constant  wheel  torque  maneuvers 
for  actual  vehicles  (GPS  Block  HR  and  Hubble  Space  Telescope).  Limitations  of  the  sub- 
optimal  control  are  explored  with  several  computer  simulations  varying  the  magnitudes  of 
the  control  torques  used. 

Chapter  5  then  introduces  gravity  gradient  and  solar  pressure  torques,  and  the  result¬ 
ing  deviations  from  the  torque  free  case  are  discussed  by  comparing  computer  simulations 
for  both  unperturbed  and  perturbed  models. 

The  overall  results  are  presented  and  discussed  in  Chapter  6,  along  with  recommen¬ 
dations  for  further  research. 
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II.  Equations  of  Motion 


In  this  chapter  we  derive  the  rotational  and  attitude  equations  of  motion  for  a  rigid 
spacecraft  TZ  with  n  embedded  axisymmetric  momentum  wheels  W„  operating  in  an  en¬ 
vironment  with  external  perturbing  torques  (solar  pressure  and  gravity  gradient).  Both 
dimensional  and  non-dimensional  forms  of  the  equations  are  developed.  The  notation  used 
is  primarily  that  of  Hughes  [11]  and  Hall  [8]. 


2. 1  Rotational  Equations  of  Motion 


Figure  2.1  depicts  a  spacecraft  composed  of  a  platform  TZ  and  momentum  wheels 


Figure  2.1  Rigid  Spacecraft  with  Multiple  Momentum  Wheels 


The  platform  is  rigid,  has  constant  mass,  is  not  necessarily  symmetric,  and  the 
body  fixed  reference  frame  Ti  with  basis  (bi  bj  baj  is  centroidal  but  not  necessarily 
principal.  The  body  fixed  frame  translates  and  rotates  relative  to  the  inertially  fixed 
reference  frame  Ei  with  basis  (ci  62  63).  The  spin  axes  of  the  momentum  wheels  are 
assumed  to  be  fixed  within  the  spacecraft.  The  wheels  are  assumed  to  be  perfectly  balanced 
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and  axisymmetric,  spinning  about  their  individual  axes  of  symmetry  ai...a„.  Before 
proceeding  with  derivation  of  the  equations  of  motion,  a  few  definitions  are  necessary. 

We  define  the  matrix  A  such  that  the  columns  are  the  column  matrices  (j  =  1 . . .  n) 
specifying  the  orientation  of  the  wheels  Wi . . .  Wn  within  the  vehicle  in 


(  : 


A  = 


ai  a.2 


\  ■■ 


:  \ 


(2.1) 


The  moments  of  inertia  matrix  for  the  entire  spacecraft  (including  wheels)  is  assumed 
constant  as  a  resvdt  of  the  constant  mass  and  rigidity  assumptions. 


Il2 

^21 
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hs 

hi 
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(2.2) 


Finally,  Is  is  a  diagonal  matrix  composed  of  the  axial  moments  of  inertia  of  the  wheels 
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We  will  now  derive  the  equations  of  motion  for  a  rigid  spacecraft  with  embedded 
momentum  wheels  for  the  general  case  in  an  environment  with  external  perturbing  torques. 

Following  Hughes  [11:159],  the  total  angular  momentum  of  a  rigid  platform  with 
embedded  momentum  wheels  can  be  expressed  in  the  form 


h  =  ic5  +  Ahs 


(2.4) 
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where  the  axial  angular  momenta  of  the  wheels  about  their  spin  axes  in  Tb  are 


h.  =  I,u>, 


(2.5) 


The  absolute  time  derivative  of  angular  momentum  (in  Ti)  can  be  written  as  the  relative 
time  derivative  of  angular  momentum  (in  Tb)  plus  the  cross  product  of  the  angular  rates 
of  Tb  relative  to  and  the  angular  momentum 


where 


^  0  -U)3  Q2  ^ 


= 


(2.6) 


U?3  0  — Ui 

-^2  wi  0  y 

In  an  environment  free  of  external  torques,  the  absolute  time  derivative  of  angular  mo¬ 


mentum  (in  J^i)  is  zero  ( 
as 


=  0),  so  that  this  equation  may  be  rearranged  and  written 


(2.7) 


L 


On  the  other  hand,  if  external  perturbing  torques  are  present,  the  absolute  time  derivative 


of  angular  momentum  (in  T,)  is  equal  to  the  sum  of  the  external  torques  ( 
that  the  equivalent  and  more  general  equation,  after  some  rearranging,  is 


L 


=  =  M-l- 


M),  so 


(2.8) 


L 


Hughes  [11]  also  shows  that  the  total  angular  momenta  of  the  wheels  is  given  by 


ho  =  I,  A^a>  -1-  h. 


(2.9) 


so  that  the  time  derivative  of  the  total  angular  momenta  of  the  wheels  can  be  expressed 
as 

k  =  ^  =  ga  (2.10) 

at 
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where  ga  are  the  torques  applied  by  the  platform  to  the  wheels  about  their  spin  axes. 


Now,  dropping  the  reference  frame  notation  (thereby  assuming  all  future  equations 
are  expressed  in  the  rotational  equations  of  motion  for  a  spacecraft  composed  of  a  rigid 
platform  and  n  embedded  rigid  momentum  wheels  Wj  {j  =  l...n)  can  be  summarized 
as 


r  dh  - _ 

h  =  —=  =  ho?  +  M 
at 

(2.11) 

i  dha 

(2.12) 

where 


h  =  id>  +  Ahj  =total  angular  momentum  of  7^  +  Wj 
ha  =  i,  +  h,  =  total  angular  momenta  of  Wj 
ii,  =  ijd;,  =  axial  angular  momenta  of  Wj  about  a^ 
M  =  external  torques  acting  on  the  vehicle 
go  =  axial  torques  applied  by  TZ  to  VV,-  about  a^- 
i  =  time 
j  =  l...n 


(2.13) 

(2.14) 

(2.15) 

(2.16) 

(2.17) 


2-4 


2.2  Non-dimensionalization 


In  order  to  more  easily  derive  the  stationary  platform  condition  and  visualize  torque 
trajectories  that  satisfy  this  condition,  it  will  be  useful  to  non-dimensionalize  the  rotational 
and  wheel  torque  differential  equations  (2.11)  and  (2.12).  This  has  the  added  advantage 
of  generalizing  the  solution  to  a  variety  of  applications.  The  following  sequence  of  steps 
is  basically  an  expansion  of  the  non-dimensionalization  developed  by  Hall  [8],  with  the 
addition  of  external  torques. 

Solving  Eq.  (2.14)  for  h,  and  substituting  into  Eq.  (2.13)  yields 


h,  =  ha  -  i,A^a) 

(2.18) 

ii  =  iw  -1-  A  ^iia  -  i,  A^w^  “  “  AIs  A^j  cb  -f  Aha 

(2.19) 

Now  define 

j  =  i  -  At  A’’ 

(2.20) 

so  that  Eq.  (2.19)  becomes 

h  =  ja>  -|-  Aha 

(2.21) 

Solving  for  results  in 

cb  =  ^h  -  Aha^ 

(2.22) 

Substituting  w  into  Eq.  (2.11)  yields 

h  =  h^  [j-'  (h-  Aha)]  -hM 

Now,  we  introduce  the  following  transformations  to  non-dimensionalize: 

11 

(2.23) 

ha 

^  "  T 

ho 

(2.24) 

II 

(2.25) 
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where  /<;  =  trace  is  a  characteristic  inertia  and  hg  is  the  magnitude  of  the  total  angular 
momentum.  Note  that  for  the  external  torque-free  case,  h  =  constant,  so  the  magnitude  at 
any  time  can  be  used.  However,  in  the  case  with  external  torques  present,  the  magnitude 
of  total  angular  momentum  can  vary.  Thus,  it  is  best  to  select  a  known  value  such  as  the 
initial  total  angular  momentum,  hg,  which  is  valid  in  both  cases. 

Since  all  moments  of  inertia  (including  J)  are  non-dimensionalized  by  dividing  by  4, 
we  can  replace  J  with  JIc  in  Eq.  (2.21)  to  obtain 


h  =  J/cd)  -I-  Aha 


(2.28) 


Dividing  both  sides  of  this  equation  by  hg  and  noting  the  relationships  given  in  Eq.  (2.23) 
and  Eq.  (2.24),  this  becomes 

X  =  -f  An  (2.29) 

hg 


Solving  this  equation  for  d>  yields 


(x  -  An) 


(2.30) 


This  is  essentially  another  derived  non-dimensionalization  transformation  relating  uj  and 
(jj  where 

(V  =  3~^  (x  -  An)  (2.31) 


We  can  thus  substitute  Eq.  (2.30)  into  Eq.  (2.11)  to  obtain 


di 


J  ^  (x  —  An) 


-l-M 


(2.32) 
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From  the  non-dimensionalization  transformation  for  i  in  Eq.  (2.27)  we  obtain 


di  =  dti^ 

w 


(2.33) 


Substituting  this  into  Eq.  (2.32)  yields 


dh  (  ho 

[z 


ho 


J  ^  (x  —  Afi) 


+  M 


(2.34) 


Multiplying  through  by  (^idhZj,  and  regrouping  slightly  results  in 


1  dh  /  h  \  1  ,  A  \  I 

ho  dt  \ho)  hi 


(2.35) 


Recognizing  from  Eqs.  (2.23)  and  (2.25)  that  h/ho  =  x,  dhjho  =  dx,  and  —  M, 

this  equation  reduces  to 

fJ'V 

^  =  x^  [J-1  (x  -  A/i)]  +  M  (2.36) 

at 


which  results  in  the  dimensionless  rotational  equations  of  motion 


X  =  x^  [J  ^  (x  -  A/t)]  +  M 


(2.37) 


or,  noting  from  Eq.  (2.31)  that  w  =  J  ^  (x  —  Afi): 


X  =  +  M 


(2.38) 


Similarly,  by  making  the  substitution  from  Eq.  (2.33)  for  di in  the  dimensional  wheel  torque 
relationship,  Eq.  (2.12)  becomes 

dhg 
dt 
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Dividing  both  sides  by  ho  and  recognizing  from  Eq.  (2.24)  that  dhafho  =  dfi  results  in 


d/r  /  ho  \  _  Qo 

^  [zj  ~  k 


(2.40) 


which  can  be  rewritten  as 


(2.41) 


Comparing  this  to  the  non-dimensionalization  transformation  relationship  for  e  given  by 
Eq.  (2.26),  we  obtain  the  non-dimensional  wheel  torque  equation 


(2.42) 
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2.3  Perturbative  Force  Modeling 

In  the  previous  sections,  we  derived  the  rotational  and  attitude  equations  of  motion 
for  a  spacecraft,  without  specifying  the  external  moments  (M).  In  real  applications,  this 
quantity  is  a  composite  of  dozens  of  external  moments  acting  on  the  vehicle.  It  is  our 
desire  to  model  the  most  significant  of  these.  In  Chapter  IV,  we  will  see  the  effects  these 
moments  have  on  the  attitude  dynamics  of  the  vehicle  during  maneuvers. 

The  external  moments  acting  on  a  spacecraft  in  orbit  about  the  Earth  vary  in  mag¬ 
nitude  (depending  on  altitude  and  spacecraft  configuration)  from  large  effects  such  as 
aerodynamic  drag  in  low  earth  orbit,  to  smaller  effects  such  as  gravitational  interactions 
between  the  vehicle  and  the  Sun,  Moon,  and  other  planets. 

Gravity  gradient  and  solar  pressure  effects  are  perhaps  the  most  significant  from 
an  attitude  control  perspective.  Deviations  from  a  pure  two-body  problem  (e.g.  Earth 
oblateness  and  third  body  effects  from  the  Sun,  Moon,  and  other  planets)  primarily  affect 
the  satellite’s  orbit.  These  effects  are  certainly  important  for  accurate  orbit  estimation. 
However,  even  though  the  inclusion  of  gravity  gradient  and  solar  pressure  torques  couples 
the  orbital  and  rotational  equations  of  motion,  the  assumption  of  an  inertially  fixed  orbit 
with  perfectly  known  orbital  elements  is  sufficient  for  preliminary  analysis. 

Accordingly,  we  make  the  assumptions  in  the  following  derivations  that  the  space¬ 
craft’s  translational  motion  is  that  of  a  point  mass  in  orbit  about  a  perfectly  spherical  Earth 
of  uniform  mass  distribution,  with  no  interaction  with  the  atmosphere  or  other  celestial 
bodies.  Orbital  parameters  for  the  GPS  Block  HR  and  Hubble  Space  Telescope  are  listed  in 
Appendices  B  and  C.  Both  vehicles  have  nearly  circular  orbits  with  altitudes  of  20, 182  km 
and  615  km,  respectively,  and  inclinations  of  54°  and  45°.  For  ease  of  simulation,  the  right 
ascension  of  the  node  (O)  was  assumed  to  be  zero  and  the  argument  of  perigee  (a;)  was 
assumed  to  be  7r/2  for  both  vehicles.  While  u)  is  not  defined  for  exactly  circular  orbits, 
the  simulation  was  left  as  general  as  possible  to  allow  application  to  elliptical  orbits.  Such 
assumptions  allowed  easy  calculation  of  initial  velocity  components  in  the  inertial  frame, 
n  and  oj  were  also  assumed  constant,  which  is  reasonable  if  the  vehicles  are  kept  on  station 
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(either  automatically  or  via  ground  commands)  or  if  simulation  times  are  limited  to  only 
a  few  orbital  periods. 


The  two-body  equations  of  motion  are  well  known  and  can  be  expressed  in  many 
forms.  However,  in  our  inertially  fixed  Earth-centered  frame  Ti,  they  can  be  expressed  as 
a  set  of  six  first  order  differential  equations  (in  matrix  form)  as 


i  =  V  (2.43) 

V  =  (2.44) 

R 

where 


=  Earth  gravitational  parameter 

R  =  components  of  position  vector  from  Earth  center  to  spacecraft  center  of  mass 
V  =  components  of  velocity  vector  of  spacecraft  center  of  mass 
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2.3.1  Gravity  Gradient  Torques.  Gravity  gradient  torques  arise  due  to  asymme¬ 
tries  in  the  vehicle’s  mass  distribution  and  the  resulting  inequity  of  gravitational  attraction 
between  the  central  attractive  body  (e.g.  Earth)  and  spacecraft  mass  elements  for  certain 
vehicle  orientations.  If  the  vehicle  were  perfectly  spherical  and  of  uniform  mass  distri¬ 
bution,  the  sum  of  gravitational  forces  on  the  differential  elements  of  mass  comprising 
the  vehicle  would  offset  each  other  resulting  in  no  net  torque  on  the  vehicle  regardless  of 
orientation.  However,  as  is  the  case  with  most  modern  spacecraft  (including  GPS  Block 
HR  and  the  Hubble  Space  Telescope),  mission  requirements  result  in  asymmetry  in  mass 
distribution. 

Following  Kaplan  [13:200]  and  Hughes  [11:235],  the  net  gravity  gradient  torque  (M^) 
acting  on  a  body  B  can  be  expressed  as 

M,  =  /  ?X 


-M® 


(i+f) 


dm 


R-l-  r 


(2.45) 


where  R  =  vector  from  attracting  body  center  of  mass  to  B  center  of  mass 
r  =  vector  from  B  center  of  mass  to  differential  mass  element  dm 

It  is  the  dependency  of  on  R  which  couples  the  orbital  equations  of  motion 
(Eq.  (2.44))  and  the  rotational  equations  of  motion  (Eq.  (2.11))  since  the  gravity  gradient 
torque  is  one  component  of  the  overall  external  moment  M.  Thus,  as  the  orbital  equations 
are  integrated,  the  resulting  position  data  are  used  to  calculate  gravitational  torques  which 
are  in  turn  used  in  the  integration  of  the  rotational  equations  of  motion. 

Summaries  of  the  mass  properties  of  the  GPS  Block  HR  and  Hubble  Space  Telescope 
vehicles  are  included  in  Appendices  B  and  C.  While  detailed  mass  data  for  every  vehicle 
subsystem  and  location  were  available  for  GPS  Block  HR  [17],  only  total  mass  was  available 
for  the  Hubble  Space  Telescope  [2:1]. 
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The  most  accurate  yet  easily  calculable  approach  to  determining  the  resultant  torques 
is  to  use  a  series  expansion  for  R  +  r  ,  recognizing  that  r  /  R  <1.  After  dropping 
squared  and  higher  terms,  and  some  manipulation  (see  Hughes  [11:235-237]),  the  expression 
for  total  gravity  gradient  torque  (in  matrix  form  in  Th)  becomes 

M,  =  ^C3^ic3  (2.46) 

R 

where  C3  =  [ci3  C23  033]^  and  Ci3,  (i  =  1 . . .  3)  are  the  direction  cosines  found  by  calculating 

b.-(i/  i). 

Analysis  of  the  above  expression  reveals  that  for  an  asymmetric  vehicle  (principle 
moments  of  inertia  unequal),  the  net  gravity  gradient  torque  on  the  vehicle  is  zero  when 
any  of  the  principal  axes  are  aligned  with  the  local  vertical  (R/  R  ).  Stability  analysis  of 
these  equilibria  shows  that  the  only  stable  orientation  for  a  vehicle  in  a  circular  orbit  is 
such  that  the  minor  axis  is  aligned  with  the  local  vertical  [11], [5], [13]. 

Figure  2.2  shows  the  relationship  between  the  body  frame  and  inertial  frame  for  the 
GPS  Block  HR  spacecraft  in  orbit.  The  relationship  for  Hubble  Space  Telescope  is  similar 
(see  Appendix  C  for  body  frame  orientation). 
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Figure  2.2  GPS  Block  HR  Gravity  Gradient  Model 


2.3.2  Solar  Pressure  Torques.  Small  torques  on  a  spacecraft  can  also  occur  as  a 
result  of  uneven  absorption  and  reflection  of  electro-magnetic  radiation  (photons)  by  the 
spacecraft  surfaces.  The  primary  source  of  high  energy  radiation  for  an  Earth-orbiting 
satellite  is  obviously  the  Sun,  and  the  effects  are  clearly  greater  for  large  asymmetric 
vehicles  than  for  smaller,  more  symmetric  ones.  Since  both  the  GPS  Block  HR  and  Hubble 
Space  Telescope  vehicles  have  large  appendages  (solar  arrays)  and  the  main  bodies  are 
asymmetric,  the  effects  of  solar  pressure  induced  torques  on  a  maneuvering  vehicle  are  of 
interest. 

As  is  more  fully  derived  in  Agrawal  [1:133-135],  the  net  solar  pressure  induced  torque 
on  a  surface  S  with  uniform  material  properties,  discounting  shadowing  effects,  is  given 
by; 

-  (n.s)  ?  X  (/>.  +  Pi)  S  -f  {2p,  +  ^)  n  (2.47) 

or,  since  Pa  +  Pd  +  =  1: 

M,  =  J  -  (n*sj  r  X  (1  -  p,)  S  -I-  2  ^Pj  -h  n  P,undA  (2.48) 
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where: 


S  =  unit  vector  from  Sun  to  differential  area  dA 
n  =  externally  oriented  unit  normal  vector  of  dA 
V  —  vector  from  vehicle  center  of  mass  to  the  center  of  pressure  of  dA 
Psun  =  solar  radiation  pressure 
Pa  =  fraction  of  radiation  absorbed  by  dA 
Pd  =  fraction  of  radiation  diffusely  reflected  by  dA 
Ps  —  fraction  of  radiation  specularly  reflected  by  dA 


It  is  important  to  restate  that  these  equations  do  not  take  into  account  shadowing 
effects  of  one  surface  on  another  (e.g.  solar  array  shadowing  main  body  or  “back”  face  of 
spacecraft  entirely  in  shadow).  For  a  three-dimensional  surface  such  as  the  main  body  of 
the  GPS  Block  HR  and  the  Hubble  Space  Telescope,  if  no  provision  is  made  to  exclude  the 
side  not  facing  the  Sun  from  the  equations,  the  resulting  torques  calculated  will  in  effect 
assume  that  the  photons  have  traveled  through  the  front  face  and  are  impinging  on  the 
rear  faces  as  well. 

Unfortunately,  the  main  body  surfaces  cannot  reasonably  be  ignored  for  either  GPS 
Block  HR  or  the  Hubble  Space  Telescope,  since  the  moment  arms  (from  vehicle  mass  center 
to  surface  center  of  pressure)  of  the  solar  pressure  torques  about  the  bi  and  b2  axes  are 
much  greater  than  those  of  the  solar  arrays,  and  the  assumption  of  symmetry  of  the  two 
solar  arrays  with  respect  to  the  bi  -  ba  plane  negates  any  torques  about  the  63  axis  from 
the  arrays. 

For  a  spacecraft  with  discrete  components,  the  solar  pressure  torque  calculation  is 
reduced  to  a  sum  of  the  contribution  from  each  of  the  component  surfaces,  recognizing 
that  any  main  body  surfaces  in  shadow  (n,-S  >  0)  contribute  nothing.  In  matrix  form, 
the  equation  is 


=  1  •- 


-(nfS)ff 


(l-p,)S  +  2(p,  +  ^ 


n. 


(2.49) 
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where: 


S  =  components  of  unit  vector  from  Sun  to  vehicle 
n,-  =  components  of  externally  oriented  unit  normal  vector  of  A, 

f,-  =  components  of  vector  from  vehicle  center  of  mass  to  the  center  of  pressure  of  Ai 

Note  in  Figure  2.3  that  the  GPS  Block  HR  vehicle  is  modeled  as  eight  planar  rect¬ 
angular  surfaces,  only  six  of  which  can  ever  be  in  shadow  (the  main  body  “sides”,  “top”, 
and  “bottom”).  For  the  Hubble  Space  Telescope  (Figure  2.4),  however,  the  two  main  body 
cylinders  are  modeled  such  that  half  the  cylinders  are  always  in  sun  unless  the  ba  axis  is 
parallel  to  S.  Due  to  our  choice  of  body  frame,  the  effective  areas  of  these  half  cylinders  are 
rectangles  with  widths  equal  to  the  cylinder  diameters  and  heights  equal  to  the  cylinder 
heights,  scaled  to  account  for  Sun  vector  angle  by  the  nf  S  term.  The  center  of  pressure 
thus  “moves”  on  the  surface  such  that  it  is  always  in  the  center  of  the  effective  area,  and 
the  normal  vector  is  defined  as  the  surface  normal  at  the  center  of  pressure  (with  only  bi 
and  b2  components). 

Differences  in  the  absorptivity  and  reflectivity  of  various  parts  of  the  spacecraft  will 
also  result  in  solar  pressure  differentials.  For  example,  it  would  be  expected  that  solar 
arrays  have  a  higher  absorptivity  than  the  main  body,  which  is  typically  covered  in  layers 
of  reflective  material.  However,  for  simplicity  of  modeling,  the  material  properties  of  both 
the  GPS  Block  HR  and  Hubble  Space  Telescope  were  assumed  constant  across  all  surfaces, 
with  spectral  reflectivity  higher  than  diffuse  reflectivity  and  absorptivity  (see  Appendices 
B  and  C  for  assumed  values). 

Due  to  the  large  distances  from  the  Sun  to  the  vehicle  for  an  Earth-orbiting  satellite 
and  the  relatively  small  size  of  the  vehicle,  it  is  reasonable  to  assume  that  the  Sun  vector 
direction  is  constant  across  all  exposed  surfaces  for  any  point  in  time.  However,  since 
the  Earth  is  moving  relative  to  the  Sun,  the  Sun  vector  does  move  with  respect  to  our 
Earth-centered  inertial  reference  frame  Since  the  Earth’s  orbit  is  nearly  circular  and 
by  definition  of  the  ecliptic  plane  has  inclination  of  zero  [3],  it  is  reasonable  to  assume  that 
the  Sun  vector  “rotates”  at  a  constant  angular  rate  in  the  ei  —  62  plane  of  our  inertial 
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reference  frame.  If  the  angle  of  the  Sun  vector  from  the  ei  axis  is  given  by  9^^  then  the 
Sun  vector  in  Ti  can  be  represented  in  matrix  form  as 


S  =  [cos  05  sin  05  0]^ 


(2.50) 


and  the  angle  05  is  governed  by 


0.1?  —  — 


27r 


365.25  days 


1.99  X  10  ”  rad/s  =  constant 


(2.51) 
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2.4  Attitude  Equations 

There  are  numerous  choices  of  parameters  available  to  express  the  attitude  of  the 
vehicle  body-fixed  frame  Eb  with  respect  to  the  inertial  frame  Each  has  its  own 
advantages  and  disadvantages  depending  on  the  application. 

Euler  angles,  usually  represented  as  three  single  axis  rotations  about  the  body  axes, 
are  perhaps  the  most  commonly  used  parameters  in  spacecraft  design,  mission  analysis, 
and  spacecraft  operations.  This  is  due  to  ease  of  visualization  and  intuitive  relationship 
between  the  body  frame  (containing  spacecraft  sensors  and  communications  hardware)  and 
the  inertial  frame  (from  which  the  sensors  obtain  data).  However,  the  well  documented 
singularities  in  the  differential  equations  relating  vehicle  angular  velocities  to  Euler  angles 
leads  to  computational  diflficulties  for  passage  through  certain  orientations  (depending  on 
the  Euler  angle  set  used). 

Attitude  equations  utilizing  quaternions,  with  a  redundant  fourth  parameter,  avoid 
these  computational  singularities.  However,  the  cost  of  numerical  stability  is  a  loss  of 
simple  physical  meaning  of  the  four  parameters  to  planners  and  analysts. 

In  all  of  the  spacecraft  command  ajid  control  and  simulation  systems  this  author  has 
dealt  with,  quaternions  are  used  for  attitude  propagation  within  the  modeling  software, 
and  Euler  angles  are  calculated  from  the  quaternions  for  display  to  the  end  user. 

The  quaternions  are  dimensionless  parameters  by  definition.  Following  Chobotov 
[5:10],  if  the  four  quaternion  parameters  are  represented  by  the  column  matrix  q  = 
(Qi  Q2  Qs  94)^)  the  attitude  equations  can  be  written  as 

0  Wa  -Cj2  Cji 
— C<^3  0  CJ\  LO2 

UJ2  — ^1  ^  ^3 

— UJl  — UJ2  — ^3  0 
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It  is  really  unnecessary  to  distinguish  between  the  dimensional  and  non-dimensional 
forms  of  this  equation,  since  if  we  substitute  in  for  dt  from  Eq.  (2.33)  and  note  the  rela¬ 
tionship  between  u)  and  a;  from  Eq.  (2.31),  Eq.  (2.52)  becomes 


^  0 

OJS 

—0^2 

dq  1 

1 

-CJa 

0 

^1 

OJ2 

dt 

Vc) 

~  2  ' 

UJ2 

-CUi 

0 

UJs 

\  -^1 

—0^3 

0 

(2.53) 


Clearly,  elimination  of  common  terms  from  both  sides  thus  yields  the  non-dimensional 
form  of  the  equation,  which  is  the  same  form  except  that  non-dimensional  angular  rates 
are  used.  Using  the  superscript  notation  to  represent  the  4x4  skew  symmetric  matrix 
of  angular  velocities  (not  to  be  confused  with  the  3  x  3  matrix  uj’^),  the  attitude  equations 
can  be  expressed  in  the  compact  matrix  form 

q=lfiXq  (2.54) 


The  quaternions  obviously  vary  with  time  if  the  spacecraft  is  rotating  with  respect 
to  the  inertial  frame.  Certainly,  the  orientation  of  one  orthogonal  reference  frame  to 
another  can  be  described  by  a  3  X  3  rotation  matrix  T  that  also  varies  with  time,  and 
whose  columns  are  the  orthonormal  basis  vectors  of  expressed  in  terms  of  the  basis 
vectors  of  Again  following  Chobotov  [5:9],  the  elements  of  this  rotation  matrix  can  be 
expressed  in  terms  of  the  quaternions  as 

Til  =  ql-ql-ql  +  ql 
Ti2  =  2(51^2  +  9394) 

Ti3  =  2  (51 53  -  5294) 

T21  =  2  (5152  -  9394) 

T22  =  -9i  +92-91  +  94 
T23  =  2  (5194  +  9293) 
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Tsi  —  2  (^1^3  +  g2Q'4) 

T32  =  2  (—91^4  +  9293) 

T22  =  -9?  -92  +  91  +  94 

Furthermore,  since  the  rotation  matrix  can  also  be  constructed  as  the  product  of 
three  single-axis  rotation  matrices  involving  Euler  angles  (choice  of  sequence  dependent  on 
application),  it  should  be  clear  that  the  Euler  angles  can  be  extracted  from  the  T  found 
using  quaternions.  Angular  rates  in  terms  of  the  Euler  angles  can  likewise  be  derived.  The 
reader  is  referred  to  Chobotov  [5:6-7]. 
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2.5  Equations  Summary 

Before  continuing,  it  will  be  useful  to  have  all  of  the  equations  of  motion  summa¬ 
rized  in  one  location  for  future  reference  (renumbered  as  well).  Note  that  only  the  non- 
dimensional  forms  of  the  rotational  and  wheel  torque  equations  are  given,  since  they  are 
used  exclusively  in  the  following  sections.  For  subsequent  developments  which  assume  no 
external  torques  (as  in  Chapter  III),  the  rotational  equations  of  motion  do  not  contain  M, 
and  thus  the  gravity  gradient  and  solar  pressure  torque  equations  used  to  calculate  the 
components  of  M  (i.e.  and  M,)  are  not  necessary.  The  following  equations  are  for 
the  general  case  with  external  torques  included  (which  will  be  used  in  Chapter  5),  and  are 
provided  in  matrix  form.  All  equations  assume  the  matrices  are  expressed  in  Tt  except 
the  orbital  position,  orbital  velocity,  and  Sun  vector  equations.  For  definitions  of  variables 
and  constants  the  reader  is  referred  to  the  previous  sections  and  to  the  List  of  Symbols. 


X  = 

x^u>  4-  M 

(2.55) 

ft  = 

e 

(2.56) 

q  = 

2  ^ 

(2.57) 

r] 

fvl 

(2.58) 

■ 

L 

v]  = 

-  [r] 

(2.59) 

R 

i  = 

^sun 

(2.60) 

M  =  ((f)(M,  +  M,) 

\ho/ 

M,  =  T^c-ics 
R 

M,  =  ^  -  (nf  S)  ff  (1  -  p,)  S  -f  2  ^p,  +  n,-  P^unA 

[S]^.  =  [cos  6,  sin0,  0]^ 


(2.61) 

(2.62) 

(2.63) 

(2.64) 
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IIL  Wheel  Torque  Control  Law  Derivation 

The  purpose  of  this  chapter  is  to  state  the  conditions  on  spacecraft  momentum 
and  wheel  momenta  necessary  for  a  stationary  platform  (body  angular  rates  very  small), 
calculate  some  wheel  torque  controls  to  accomplish  rest-to-rest  and  stationary  platform 
maneuvers  and  discuss  limitations  on  such  maneuvers  due  to  physical  limitations  of  the 
wheels. 

Examples  of  these  maneuvers  for  specific  spacecraft  will  be  performed  in  Chapter 
4  (without  external  perturbing  torques)  and  Chapter  5  (with  gravity  gradient  and  solar 
pressure  torques). 

3,1  Stationary  Platform  Condition 

Most  spacecraft  in  operation  today  are  used  for  communications  or  observation  and 
thus  have  mission  constraints  that  require  them  to  remain  Earth  pointed  (e.g.  GPS)  or 
essentially  fixed  in  inertial  space  (e.g.  Hubble  Space  Telescope)  for  long  periods  of  time. 
This  necessitates  the  ability  to  maintain  the  spacecraft  attitude  and  counteract  external 
perturbations  such  as  gravity  gradient  and  solar  pressure  torques.  Also,  especially  in  the 
case  of  the  Hubble  Space  Telescope,  it  is  often  critical  to  be  able  to  maneuver  the  vehicle 
from  one  stationary  orientation  to  another  with  minimal  induced  motion  on  the  vehicle 
(which  might  damage  components  such  as  solar  arrays,  antennas,  or  sensors). 

It  is  immediately  apparent  from  inspection  of  Eq.  (2.31)  that  the  vehicle  will  be 
stationary  if 

10  =  J-^(x-A/x)  =  0  (3.1) 

which  leads  to  the  stationary  platform  condition 

X  =  Afi  (3-2) 


For  vehicles  with  3  non-coplanar  wheels,  this  can  also  be  expressed  as 

/i  =  A~^x 


(3.3) 
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The  physical  interpretation  of  Eq.  (3.2)  is  quite  simple;  for  the  platform  to  be  sta¬ 
tionary,  the  direction  and  magnitude  of  the  spacecraft  total  momentum  vector  (x)  must 
be  equal  to  the  direction  and  magnitude  of  the  total  wheel  momentum  vector  ( A/x)  in  the 
body  frame  for  all  time. 

An  interesting  way  of  describing  and  visualizing  this  stationary  platform  condition 
was  developed  by  Hall  [8]  by  noting  that  since  x  is  normalized,  conservation  of  angular 
momentum  requires  x^x  =  1,  so  that  substituting  from  above  ultimately  yields 

fi^A^AfjL  =  1  (3.4) 

which  for  a  vehicle  employing  n  momentum  wheels,  is  readily  recognized  as  a  hyper-ellipsoid 
of  dimension  n  in  the  space  spanned  by  (fii  0  •  •  -  0),  (0  ^2  0  ■  ■  -  0),  (0  ■  •  -  0  /x„),  hereafter 
referred  to  as  /x  space  (or  Thus,  for  a  given  value  of  total  angular  momentum  (ho), 
the  platform  TZ  will  remain  fixed  with  respect  to  inertial  space  {Ti)  as  long  as  the  wheel 
momenta  fi  =  {pi  P2  •  •  •  Pn)  lie  on  the  surface  of  this  wheel  momenta  ellipsoid  and  x  =  A/x. 

Also,  as  is  common  in  dynamics  texts  and  attitude  dynamics  papers,  the  spacecraft 
total  angular  momentum  itself  can  be  represented  by  a  sphere,  in  our  case  in  x  space  (see 
Fig.  3.1).  Typically,  this  momentum  sphere  is  used  to  visualize  the  stability  of  spacecraft 
motion  about  the  angular  momentum  vector.  However,  in  our  applications,  we  will  use  this 
momentum  sphere  only  to  describe  the  spacecraft  total  momentum  during  reorientation 
maneuvers.  It  will  prove  particularly  useful  in  showing  how  violations  of  assumptions  on 
control  torque  magnitudes  results  in  state  errors. 

While  only  three  wheels  are  necessary  for  three  degree  of  freedom  attitude  control, 
most  vehicles  employing  them  have  a  set  of  four,  aligned  in  a  “pyramid”  fashion  so  that  no 
combination  of  three  of  the  wheels  are  co-planar  (see  Appendix  C  for  GPS  Blk  HR  wheel 
orientations).  Having  a  fourth  wheel  not  only  provides  redundancy  in  the  event  of  failure 
or  degradation  of  performance  of  a  single  wheel,  but  also  allows  more  total  momentum 
storage  by  the  wheels.  Thus,  while  in  some  spacecraft  only  three  of  the  four  wheels  are  used 
simultaneously  (e.g.  Hubble  Space  Telescope),  in  others,  the  fourth  wheel  is  routinely  used 
to  increase  control  authority  (e.g.  GPS  Block  HR)  [20],[2],[7].  For  ease  of  visualization. 
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Figure  3.1  Spacecraft  Total  Momentum  Sphere 

we  will  focus  primarily  on  reorientation  maneuvers  using  sets  of  three  wheels.  This  allows 
representation  of  the  stationary  platform  condition  by  a  three-dimensional  wheel  momenta 
ellipsoid.  However,  in  later  sections  we  will  briefly  discuss  use  of  the  fourth  wheel  in  cases 
where  the  primary  three  wheels  cannot  generate  enough  angular  momentum  to  maintain 
a  stationary  platform. 

Figures  3.2  through  3.9  are  examples  of  wheel  momenta  ellipsoids  for  vehicles  using 
three  momentum  wheels  for  attitude  control,  for  varying  orientations  of  the  wheels  within 
the  body  frame.  Notice  that  for  the  case  of  orthogonal  wheels  (Figs.  3.2  and  3.3),  the 
ellipsoid  reduces  to  a  sphere,  and  that  as  the  wheels  approach  co-planarity  (Figs.  3.6 
and  3.7),  the  ellipsoid  elongates.  If  the  three  wheels  are  co-planar  (but  none  co-linear), 
the  principal  axis  of  the  ellipse  perpendicular  to  the  plane  becomes  infinitely  large.  This 
reflects  the  fact  that  the  wheels  can  no  longer  “absorb”  the  spacecraft’s  angular  momentum 
(thus  achieving  stationary  platform  condition),  no  matter  how  small,  if  the  direction  of  the 
angular  momentum  is  not  within  the  plane  formed  by  the  three  wheels.  As  the  three  wheels 
become  co-linear  (see  Fig.  3.9),  two  of  the  three  principal  axes  become  infinitely  large, 
suggesting  that  the  wheels  can  maintain  a  stationary  platform  only  if  the  total  angular 
momentum  vector  is  aligned  with  the  wheel  spin  axes. 
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It  is  important  to  note  that  the  discussion  in  this  chapter  will  be  limited  by  the 
assumption  that  the  spacecraft  total  angular  momentum  is  constant  (h„  =  const),  so  that 
the  momentum  sphere  and  wheel  momenta  ellipsoid  are  of  constant  “size”.  In  other  words, 
since  these  two  graphical  tools  are  constructed  in  non-dimensional  variables  (x  and  /i),  if 
the  magnitude  of  the  total  angular  momentum  changes  with  time,  the  wheel  momenta 
ellipsoid  will  “expand”  or  “contract”  relative  to  the  nominal  ellipsoid  calculated  based 
on  ho.  This  will  be  discussed  in  more  depth  in  Chapter  5  when  perturbing  torques  are 
introduced  into  the  simulations. 


Figure  3.3  Wheel  Momenta  Ellipsoid 
Figure  3.2  Orthogonal  Wheels  (Orthogonal  Wheels) 
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Figure  3.5  Wheel  Momenta  Ellipsoid 
Figure  3.4  Pyramidal  Wheels  (Pyramidal  Wheels) 


Figure  3.7  Wheel  Momentum  Ellipsoid 

Figure  3.6  Nearly  Co-planar  Wheels  (Nearly  Co-planar  Wheels) 


Figure  3.9  Wheel  Momentum  Ellipsoid 

Figure  3.8  Nearly  Co-linear  Wheels  (Nearly  Co-linear  Wheels) 


Another  important  relationship  to  understand  is  that  between  the  direction  of  the 
total  angular  momentum  vector  (x),  the  direction  of  the  wheel  momenta  vector  (Aft), 
and  the  attitude  of  the  spacecraft  with  respect  to  Let  us  assume  that  the  wheels  can 
absorb  all  the  angular  momentum  of  the  spacecraft,  so  that  the  the  stationary  platform 
condition  x  =  Afi  (which  is  expressed  in  iFi,)  holds.  Also  assume  that  we  know  (i.e.  can 
accurately  estimate)  the  direction  of  the  spacecraft  angular  momentum  vector  in  inertial 
space  (x,),  the  wheel  angular  rates  (u>,),  inertias  (Ij),  and  the  orientation  of  the  wheels 
with  respect  to  the  body  frame  (A).  Finally,  let  the  attitude  of  the  spacecraft  body  frame 
with  respect  to  inertial  be  given  by  the  rotation  matrix  T.  We  can  thus  calculate 
the  spacecraft  total  momentum  in  the  body  frame  as 

X  =  Txi  (3.5) 


Substituting  Afi  in  for  x  and  solving  for  /x  we  obtain 

fi  =  A“^Tx, 


(3.6) 
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which  is  strictly  only  valid  for  a  vehicle  with  three  non-coplanar  wheels  (A  square  and 
invertible).  Since  x,  is  constant  in  the  absence  of  external  torques,  the  vehicle  can  thns  be 
maintained  stationary  at  any  attitude  by  varying  the  wheel  momenta  /t. 

A  useful  extension  of  this  result  is  that  given  a  desired  attitude  for  the  vehicle,  known 
wheel  inertias  and  orientations  in  and  a  known  constant  total  momentum  in  we 
can  calculate  the  required  wheel  angular  velocities  to  maintain  a  stationary  platform  as 
follows: 

From  Eqs.  (2.14)  and  (2.15),  we  know 

h,  =  (3.7) 

Noting  that  for  a  stationary  platform,  d>  =  0,  we  can  solve  for  ci>s  to  obtain 

W,  =  ir'h.  (3.8) 

Recognizing  from  Eq.  (2.24)  that  ha  =  we  can  substitute  into  the  above  to  obtain 

(3-9) 

Substituting  from  Eq.(3.3)  for  yields 

w,  =  ha  (ae)  X  (3.10) 

Finally,  substituting  from  Eq.  (3.5)  for  x  results  in 

ci,  =  ha  (Ai,)“'Txi  (3.11) 

A  final  useful  extension  of  Eq.  (3.5)  is  calculation  of  final  wheel  momenta  for  a 
maneuver.  Given  rotation  matrices  relating  the  body  frame  to  inertial  at  an  initial  and 
final  time  (T^  and  T/),  initial  wheel  momenta  and  using  the  fact  that  Xj  is  constant, 
it  is  easily  shown  that 

ti^  =  A-^TjT'lAtia  (3.12) 
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As  a  simple  example,  consider  a  vehicle  with  three  orthogonal  wheels  (A  =  A~^  =  I) 
and  total  angular  momentum  vector  pointing  in  the  ei  direction  (x,-  =  (1  0  0)^).  We 

assume  the  vehicle  is  stationary  in  . 

If  we  wish  to  have  the  spacecraft  bi  axis  aligned  with  the  angular  momentum  vector 
(and  thus  parallel  to  Ci),  then  =  J,  so  by  Eq.  (3.6),  /a  =  (1  0  0)^.  In  other  words, 
to  maintain  this  attitude,  only  >Vi  is  rotating,  such  that  it  contains  all  of  the  angular 
momentum  of  the  vehicle. 

To  demonstrate  that  this  attitude  is  not  unique  and  that  any  is  achievable  (given  our 
previous  assumptions),  let  us  now  assume  that  the  vehicle  has  been  reoriented  (without 
external  forces)  by  a  negative  45°  rotation  about  the  ba  axis  and  is  again  stationary.  The 
total  momentum  has  not  changed  with  respect  to  .Fj.  Using  Eq.  (3.12),  we  find  that  now 
/i  =  (I/V2  l/-\/2  0)^,  meaning  that  both  >Vi  and  W2  are  in  use  and  each  contribute 

half  the  total  angular  momentum  of  the  vehicle.  The  angular  rates  of  the  wheels  necessary 
for  given  wheel  inertias  and  vehicle  total  angular  momentum  could  be  calculated  using 
Eq.  (3.10). 

One  final  but  important  point  which  was  identified  by  Hall  [8]  is  the  following:  Al¬ 
though  the  relationship  between  the  platform  and  momentum  wheel  angular  momentum 
vectors  direction  and  magnitude  is  uniquely  determined  by  the  stationary  platform  condi¬ 
tion  X  =  A/i,  the  attitude  of  the  spacecraft  with  respect  to  inertial  space  is  determined 
only  within  a  rotation  about  the  angular  momentum  vector.  Thus,  for  all  the  maneuvers 
we  will  perform  in  the  following  sections,  even  though  we  can  control  the  direction  of  the 
spacecraft  angular  momentum  vector  within  the  body,  we  cannot  control  the  body  rotation 
with  respect  to  that  vector.  This  will  be  particularly  important  when  analyzing  final  state 
errors  at  the  completion  of  maneuvers. 

We  will  now  address  some  additional  issues  regarding  use  of  the  wheel  momentum 
ellipsoid  to  specify  the  stationary  platform  condition,  derive  two  classes  of  wheel  torque 
control  laws  for  rest-to-rest  attitude  maneuvers,  and  run  several  representative  simulations 
for  actual  spacecraft  to  characterize  the  utility  and  behavior  of  these  control  laws. 
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3.2  Unachievable  Regions  on  the  Wheel  Momenta  Ellipsoid 

Since  the  wheel  momenta  ellipsoid  is  normalized  (scaled  by  hg),  interpretation  of 
points  on  the  ellipsoid  can  be  misleading  when  used  to  assess  momentum  transfer  in  actual 
spacecraft.  There  may  in  fact  be  regions  on  the  wheel  momenta  ellipsoid  that  are  un¬ 
achievable  because  they  require  unrealistically  high  angular  momentum  for  the  individual 
wheels.  Momentum  wheels  obviously  have  fixed  axial  inertias  and  maximum  achievable 
angular  rates,  thereby  limiting  the  total  momentum  an  individual  wheel  can  generate,  as 
well  as  the  total  momentum  that  can  be  generated  by  the  set  of  three  wheels  in  certain 
directions  within  the  body  frame  Ej,. 

Consider  a  spacecraft  using  three  orthogonal  momentum  wheels,  with  axial  inertias 
of  /j  =  0.1  kg  and  maximum  achievable  angular  rates  of  =  ±100  radjs.  The 

maximum  momentum  that  a  single  wheel  can  store  (in  its  axial  direction)  is  thus  h^ax  = 
h^smax  =  ±10  %  m?  fs.  The  wheel  momenta  ellipsoid  representing  the  stationary  platform 
condition  for  this  vehicle  is  identical  to  Fig.  3.3.  If  we  consider  two  wheels  with  spin 
axes  aligned  with  the  bi  and  b2  body  axes,  the  total  momentum  generated  by  these 
two  wheels  at  maximum  spin  rate  is  hmax  ^  ±14.14  kg  m^/s  in  the  ±(l/-\/2)bi  ± 
(l/v^)b2  direction.  Obviously,  extending  this  argument  to  three  wheels,  approximately 
±17.32  kg  m? /s  can  be  generated  in  ±(l/-\/3)bi  ±  (l/>/3)b2  ±  (l/\/3)b3  direction. 

First,  let  us  suppose  that  we  wish  to  maintain  the  platform  in  an  inertially  fixed 
orientation  and  the  spacecraft  total  angular  momentum  magnitude  is  ho  =  12  kg  m^fs, 
perhaps  the  result  of  buildup  of  angular  momentum  due  to  solar  pressure  and  gravity 
gradient  torques  over  time  (the  motivation  for  Chapter  5).  Given  the  calculations  in  the 
last  paragraph,  it  is  apparent  that  a  single  wheel  alone  cannot  produce  enough  angular 
momentum  to  maintain  a  stationary  spacecraft.  However,  two  wheels  working  together 
can.  Thus,  the  only  unachievable  regions  on  the  wheel  momenta  ellipsoid  are  at  points 
“outside”  the  six  planes  which  intersect  the  ellipsoid  at  the  point  on  the  p,  axes  where  the 
wheel  momentum  maxima  are  reached  (see  Fig.  3.10).  There  are  points  in  the  pi  —  p2i 
fJ'i  ~  fJ-si  ^2  ~  /^3  planes  that  can  be  attained,  and  thus  there  are  certain  attitudes  at 
which  the  spacecraft  can  be  held  stationary. 
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Figure  3.10  Example  1  -  Achievable  Regions  on  Wheel  Momenta  Ellipsoid 

In  the  second  example,  assume  that  the  total  angular  momentum  magnitude  is  in¬ 
creased  to  ho  =  15  kg  m?  js.  Now,  not  only  is  a  single  wheel  unable  to  maintain  a  stationary 
platform,  but  any  two  wheels  working  together  are  also  insufficient.  Certainly,  three  wheels 
at  maximum  angular  rates  can  generate  enough  momentum.  In  terms  of  the  wheel  mo¬ 
menta  ellipsoid,  this  means  that  no  points  on  the  ellipsoid  in  the  Pi—p2,  Pi  —Ps-,  or  P2  —  Ps 
planes  can  be  reached  given  the  wheel  physical  limitations,  and  in  fact  only  certain  regions 
off  these  planes  are  achievable.  The  achievable  regions  on  the  wheel  momenta  ellipsoid 
are  those  regions  of  the  ellipsoid  “inside”  the  six  planes  that  define  the  maximum  wheel 
momenta  (see  Fig.  3.11). 

As  a  final  extreme  example,  consider  the  case  where  ho  =  20  kg  m^fs.  It  should  be 
obvious  that  since  even  three  wheels  working  together  at  maximum  angular  rates  cannot 
generate  more  than  17.32  kg  m^/s,  there  is  no  way  that  the  vehicle  can  be  kept  stationary 
using  the  momentum  wheels  alone  at  any  attitude.  Graphically,  the  six  planes  intersecting 
with  the  ellipsoid  mentioned  previously  have  moved  inward  so  that  the  ellipsoid  surface  is 
entirely  outside  the  cube  enclosed  by  the  planes. 
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Figure  3.11  Example  2  -  Achievable  Regions  on  Wheel  Momenta  Ellipsoid 

It  is  in  such  a  case  as  this  that  the  addition  of  the  fourth  momentum  wheel  proves 
useful.  Since  the  redundant  wheel  is  usually  oriented  such  that  it  is  non-coplanar  with  any 
other  two  wheels,  it  can  provide  a  component  of  angular  momentum  in  the  direction  of  each 
of  the  other  three  wheels,  thereby  effectively  increasing  their  momentum  storage  capacity. 
This  in  turn  allows  the  wheels  to  maintain  a  stationary  platform  at  higher  total  angular 
momentum  magnitudes.  In  this  example,  if  the  fourth  wheel  spin  axis  were  described 
in  the  body  frame  as  a4  =  (l/\/3  l/^/3  l/-\/3)^,  the  maximum  angular  momentum 

capacity  in  the  direction  of  each  of  the  other  three  primary  wheels  is  increased  by  a  factor 
of  1/y/S  w  0.58  to  approximately  15.8  kg  m?  js.  Furthermore,  in  the  direction  of  the 
fourth  wheel,  the  total  momentum  capacity  of  the  four  wheels  at  maximum  angular  rates 
is  increased  to  approximately  17.32+10.00  =  27.32  kg  m?/s.  Thus,  in  our  extreme  example 
above,  the  addition  of  the  fourth  wheel  would  result  in  regions  of  the  momenta  ellipsoid 
that  are  achievable  (similar  to  Fig  3.11),  and  thus  attitudes  for  which  the  vehicle  can  be 
held  stationary. 


3-11 


The  implications  of  the  fourth  wheel  in  maneuvers  are  apparent;  if  we  wish  to  reorient 
the  spacecraft  from  one  stationary  condition  to  another,  the  wheel  torque  control  must 
be  such  that  the  wheel  momenta  are  within  the  achievable  regions  defined  by  the  cube 
discussed  above.  Also,  if  we  choose  a  torque  control  law  that  runs  tangent  to  the  surface  of 
the  wheel  momenta  ellipsoid  (which  is  the  point  of  Section  3.5),  we  obviously  must  choose 
that  path  such  that  it  avoids  the  unachievable  regions  on  the  surface. 

Although  the  above  examples  were  for  sets  of  orthogonal  wheels,  it  should  be  apparent 
that  the  same  arguments  hold  for  non-orthogonal  wheels.  The  only  difference  is  that  since 
the  wheel  momenta  ellipsoids  are  not  symmetric  with  respect  to  the  ^  axes,  the  intersection 
of  the  six  planes  with  the  ellipsoid  produces  slightly  different  results  (see  Fig.  3.12). 

mu3 
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Figure  3.12  Achievable  Regions  on  Wheel  Momenta  Ellipsoid  (Non-orthogonal  Wheels) 
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Now  that  we  understand  the  limits  of  the  stationajy  platform  condition  and  the  use 
of  the  wheel  momenta  ellipsoids  to  visualize  it,  let  us  develop  some  torque  control  laws 
to  allow  us  to  reorient  the  vehicle  from  one  stationary  platform  state  to  another  (rest-to- 
rest).  In  Section  3.4  we  address  a  simple  constant  wheel  torque  case  to  use  as  a  baseline 
for  comparison  with  a  sub-optimal  control  law  developed  in  Section  3.5. 

Before  continuing,  however,  it  is  important  to  step  back  briefly  and  review  some 
of  the  work  accomplished  by  Hall  [8]  in  terms  of  stability  of  motion  and  allowable  wheel 
torques. 
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3.3  Allowable  Torque  Magnitudes 

Hall  [8]  shows  that  for  our  idealized  rigid  spacecraft  with  embedded  momentum 
wheels,  the  system  Hamiltonian,  using  non-dimensional  variables,  is  given  by 

H  —  -x^J“^x  -  -|-  f{C)  (3.13) 

2 

where  /(C)  is  an  arbitrary  function  of  the  first  integral  C  =  x^x/2  =  1/2  and  does  not 
alfect  the  resulting  equations  of  motion. 

Furthermore,  Hall  shows  that  the  Hamiltonian  satisfies 

it  =  =  -e^A^J-'x  (3.14) 

Ofl 

He  then  continues  by  analyzing  equilibrium  motions  of  the  integrable  zero  torque 
case  (e  =  0),  plotting  equilibrium  surfaces  in  “/tH  space”  (which  is  4-dimensional  for  a 
vehicle  with  three  wheels),  and  using  the  method  of  averaging  to  show  that  for  the  case 
of  small  constant  torques  {e  =  ea,  e  <  1,  cr  =  const)  the  equilibrium  surfaces  contain 
exact  solutions  to  the  averaged  equations.  Extending  this  approach  to  that  of  small  but 
not  necessarily  constant  torques  (such  as  the  sub-optimal  control  laws  in  Section  3.5),  he 
states  [9:10] 

. .  .trajectories  which  start  on  or  near  an  equilibrium  surface  in  [fiH]  space  will 
remain  near  the  equilibrium  surface  as  long  as  no  resonance  zone  (instantaneous 
separatrix)  is  traversed.  If  the  initial  condition  is  also  near  a  stationary  platform 
equilibrium  (a?  =  0)  and  the  torques  are  chosen  so  that  the  trajectory  in  [/i] 
space  satisfies  the  stationary  platform  condition  . . . ,  then  the  angular  velocity 
of  the  platform  will  remain  small  throughout  the  maneuver. 

Thus,  in  our  applications,  for  the  stability  analysis  to  be  strictly  valid,  all  wheel 
torque  trajectories  must  start  on  an  equilibrium  surface  in  fiH  space,  must  remain  small, 
and  should  not  cross  a  resonance  zone  (which  might  result  in  capture  into  unsteady  motion). 

Since  the  primary  purpose  of  this  report  is  to  assess  the  utility  and  limitations  of 
a  sub-optimal  wheel  torque  control  law  for  real  vehicles,  and  because  the  visualization 
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of  equilibrium  surfaces  in  fiH  space  for  a  three  wheel  vehicle  is  difficult,  stability  of  the 
various  trajectories  will  in  general  not  be  addressed. 

Stability  issues  aside,  as  we  shall  see  in  more  detail  in  the  examples  in  Chapter  4, 
violations  of  the  small  torque  assumption  also  result  in  errors  in  the  final  vehicle  states. 
Specifically,  any  movement  along  the  stationary  platform  wheel  momenta  ellipsoid  must 
violate  the  stationary  platform  condition  in  Eq.  (3.2)  because  for  there  to  be  motion  of 
the  vehicle,  no  matter  how  small,  a?  cannot  be  zero.  Thus,  the  angular  momentum  vector 
X  will  move  with  respect  to  the  platform  during  maneuvers,  and  we  shall  see  that  the 
character  of  this  motion  (large  or  small  platform  angular  velocities)  is  dependent  on  the 
choice  of  the  torque  scaling  parameter  e. 

Two  other  basic  assumptions  for  this  analysis  are  that  at  the  initial  time  the  platform 
is  stationary  and  at  the  completion  of  maneuvers,  all  torques  are  turned  off  such  that  the 
momentum  wheels  remain  at  the  final  angular  velocities  (no  energy  loss  in  wheels).  In  real 
applications,  of  course,  a  controller  using  state  feedback  may  be  necessary  to  stabilize  the 
platform  prior  to  beginning  such  a  maneuver  as  well  as  to  reduce  final  state  errors,  and 
there  may  be  bearing  or  other  losses  in  the  wheels  requiring  small  torques  to  maintain 
constant  angular  velocities  after  the  maneuver. 
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S.4  Direct  Rest-to-rest  Trajectories 

As  mentioned  previously,  there  are  certainly  applications  for  “large  angle”  spacecraft 
maneuvers  that  start  and  end  with  a  stationary  platform  (a>  «  0  at  to  and  tf).  This 
general  class  of  maneuver  is  commonly  referred  to  as  rest-to-rest,  and  includes  both  the 
arbitrary  intermediate  state  sub-class  ((i  7^  0,  t  <  t/)  and  the  stationary  platform 
maneuver  sub-class  for  which  platform  angular  rates  are  small  throughout  the  maneuver. 
In  this  section,  we  will  focus  on  a  special  case  of  the  former  sub-class,  known  as  the  direct 
rest-to-rest  maneuver.  The  latter  sub-class  will  be  discussed  in  the  next  section  of  this 
chapter. 

The  direct  rest-to-rest  maneuver  (hereafter  referred  to  as  the  direct  trajectory)  is 
depicted  in  Fig.  3.13  in  terms  of  the  wheel  momenta  ellipsoid. 


Figure  3.13  Direct  Trajectory  on  Wheel  Momenta  Ellipsoid 

It  is  apparent  that  while  both  the  initial  and  final  wheel  momenta  states  satisfy 
the  stationary  platform  criteria  (are  on  the  ellipsoid  surface),  the  intermediate  states  do 
not,  and  we  thus  do  not  expect  to  have  small  platform  angular  velocities  throughout  the 
maneuver.  As  we  shall  see  in  later  sections,  violation  of  the  stationary  platform  condition 
during  these  maneuvers  also  results  in  final  state  errors  (in  x),  so  that  even  though  the 
total  angular  momentum  is  constant  (no  external  torques),  we  can  strictly  only  use  the 
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method  outlined  in  section  3.1  to  calculate  the  attitude  of  the  vehicle  at  the  initial  and 
final  times  if  the  vehicle  is  totally  stationary. 

The  wheel  torque  control  law  for  this  direct  trajectory  from  to  is 


(3.15) 

where  e  is  the  small  parameter  discussed  in  section  3.3.  Note  that  using  this  approach,  e 
is  a  positive  constant. 

Assuming  to  =  0,  the  final  (non-dimensional)  time  for  such  a  maneuver  given  initial 
and  final  rotor  momenta  is  found  to  be 


Using  Eq.  (2.27),  the  dimensional  time  is  easily  calculated  and  is  given  by 


^  U,  J  € 


(3.17) 


By  solving  Eq.  (3.16)  for  e  and  substituting  in  for  if  from  Eq.  (3.17),  we  can  find  the  value 
of  the  small  parameter  €  required  to  accomplish  the  maneuver  in  time  if 


(k\ 

\ 


(3.18) 


Note  that  e  is  directly  proportional  to  and  while  it  is  inversely  pro¬ 

portional  to  ho  and  if.  Thus,  we  would  expect  larger  required  torques  for  larger  vehicle 
inertias  and  longer  torque  trajectories,  while  increasing  wheel  momenta  and  final  times 
will  result  in  smaller  required  torques. 

A  plot  of  e  vs.  tf  for  the  Hubble  Space  Telescope,  with  4  «  1.8  X  10®  kg  rn?,  using 
wheels  VWi  -  W3  with  =  [1  0  0],  ftf  =  [0  1  0]  for  both  ho  =  100  and  200  kg  m^/s,  is 
shown  in  Fig.  3.14. 
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Final  time  (sec)  x  1 0'* 

Figure  3.14  Direct  Trajectory  -  Small  Torque  Parameter  e  vs.  ij  for  varying  ho 

As  we  would  expect,  reducing  Ijho  by  1/2,  either  by  doubling  ho  or  halving  /<,  results 
in  required  torques  for  a  given  ij  being  reduced  by  a  factor  of  2  as  well.  Recalling  from 
our  earlier  discussion  about  allowable  torques  and  the  resulting  state  errors  during  and  at 
the  end  of  maneuvers,  we  see  that  there  is  a  significant  trade-off  between  final  time  and 
wheel  torques.  This  trade-off  will  be  discussed  in  more  detail  for  the  individual  vehicles  of 
interest  in  Chapter  4. 
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3.5  Sub-optimal  Trajectories  on  the  Wheel  Momenta  Ellipsoid 

While  the  direct  trajectories  discussed  in  the  previous  section  are  certainly  easy  to 
calculate  and  simple  to  implement,  it  will  become  apparent  from  examples  in  Chapter 
4  that  the  relatively  large  body  angular  rates  that  result  during  the  maneuver  for  large 
wheel  torques  could  prove  harmful  to  flexible  and/or  sensitive  components  on  board.  We 
therefore  wish  to  derive  a  wheel  torque  control  law  that  not  only  reorients  the  vehicle  from 
one  stationary  orientation  to  another,  but  that  maintains  small  angular  velocities  of  the 
platform  throughout  the  maneuver.  One  could  consider  such  a  maneuver  to  be  “smooth”. 

Large  angle  reorientations  for  spacecraft  using  momentum  wheels  are  often  accom¬ 
plished  via  “Euler-axis”  or  “eigen-axis”  rotations,  which  are  finite  angle  rotations  about  a 
body  fixed  axis  of  rotation  [11:10-11].  While  simple  in  concept,  calculations  to  determine 
such  maneuvers  are  somewhat  computationally  intense,  requiring  the  solution  to  an  eigen¬ 
value  problem.  Today,  such  calculations  are  typically  performed  by  high  speed  ground 
segment  computers  based  on  down-linked  telemetry,  then  sent  to  the  relatively  “dumb” 
satellite  via  the  command  up-link  for  immediate  or  delayed  execution.  However,  the  logical 
and  inevitable  progression  toward  greater  spacecraft  autonomy  (driven  by  lowered  opera¬ 
tions  costs)  dictates  on-board  attitude  determination  as  well  as  large  angle  maneuvering 
capability.  Also,  although  these  maneuvers  can  be  accomplished  rather  quickly,  they  can 
also  result  in  undesirably  large  body  angular  velocities. 

The  sub-optimal  control  law  discussed  herein  is  a  “smooth”  maneuver  in  the  sense 
that  it  reduces  unwanted  body  motions,  and  requires  only  matrix  multiplications  and 
transpose  operations  to  perform  (with  the  exception  of  final  time  calculations  which  need 
not  be  done  on-board). 

Recalling  the  discussion  in  Section  3.2,  Hall  [8]  showed  that  for  a  spacecraft  using 
two  momentum  wheels,  if  a  wheel  torque  trajectory  was  calculated  that  stayed  on  a  stable 
surface  in  the  p,H  space,  used  small  torques  (e  <  1),  and  was  tangent  to  the  stationary 
platform  wheel  momenta  ellipse,  a  “smooth”  reorientation  could  be  performed  with  only 
small  induced  body  angular  velocities. 
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Hall  derived  this  torque  control  law  by  differentiating  the  wheel  momenta  ellipsoid 
Eq.  (3.4)  with  respect  to  time,  yielding 


(3.19) 


Since  A^A  is  positive  definite  and  we  are  not  interested  in  the  trivial  solution  /x  =  0,  this 
requires  that  lie  in  the  left  nullspace  of  In  the  two  wheel  case,  it  is  thus  evident 

that  the  wheel  torque  law  is  given  by 
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which  is  easily  seen  as  a  slowly  moving  trajectory  along  the  tangent  to  the  wheel  momenta 
ellipse  in  two  dimensions.  Note  as  well  that  afai  =  1  and  a^as  =  1. 

To  extend  this  simple  result  to  three  wheels,  for  which  the  wheel  momenta  ellipsoid 
is  three  dimensional,  we  define  a  new  orthonormal  frame  with  the  same  origin  as 
but  rotated  such  that  the  i/i  axis  is  co-linear  with  and  both  and  /Xy  lie  in  the  i/i  —  P2 
plane.  The  intersection  of  the  Vi  —  z/2  plane  with  the  surface  of  the  wheel  momenta  ellipsoid 
is  thus  a  wheel  torque  trajectory  that  will  transfer  momentum  from  /x^  to  (See  Fig. 
3.15).  Note  that  strictly  there  are  two  trajectories  between  these  points,  depending  on 
the  sign  of  e  chosen  (the  “long”  and  “short”  trajectories).  We  will  be  primarily  concerned 
with  the  “short”  trajectory  since  we  expect  this  to  require  less  maneuvering  time.  For  the 
sub-optimal  control  law  derived  below,  this  corresponds  to  negative  e. 

For  a  vehicle  with  orthogonal  wheels  (wheel  momenta  ellipsoid  is  a  sphere),  this 
is  seen  to  be  the  “great  circle”  trajectory,  and  is  the  minimum  distance  (in  /x  space) 
between  /x^  and  /x^  on  the  sphere.  For  vehicles  whose  wheels  are  more  co-planar  or  co- 
linear,  this  path  obviously  is  no  longer  minimum  distance,  but  is  a  close  approximation. 
In  an  extreme  case  of  the  ellipsoid  becoming  quite  long  and  thin  (see  Fig.  3.6),  the 
minimum  distance  trajectory  would  best  be  approximated  by  the  ellipse  formed  by  taking 
a  cross-section  through  /x^,  /x^-,  and  their  reflections  across  the  longest  principal  axis  of 
the  ellipsoid.  However,  since  spacecraft  typically  have  pyramidal  wheel  arrangements. 
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A 


Figure  3.15  Sub-optimal  Trajectory  on  Wheel  Momenta  Ellipsoid 


the  wheel  momenta  ellipsoids  are  more  spherical  than  cylindrical,  and  the  simpler  result 
obtained  below  with  this  assumption  should  be  adequate. 

Finding  the  rotation  matrix  relating  Tj,  to  is  straightforward.  The  basis 
vectors  for  Tj,  expressed  in  are 
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This  orthonormal  matrix  allows  us  to  transform  a  vector  u  expressed  in  to  its  equivalent 
form  ft  in  through  the  relationship 


ft  = 


(3.25) 
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Substituting  from  this  equation  for  /t  in  Eq.  (3.4)  yields  the  stationary  platform  wheel 
momentum  ellipsoid  expressed  in 


(3.26) 
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Looking  at  the  trajectory  with  respect  to  the  new  frame  Tv,  we  note  that  the  i/3  component 
oi  V  —  {vi  V2  U3)  along  the  trajectory  is  zero,  so  that  Eq.  (3.26)  in  expanded  form  now 
reduces  to 

(  rr  rr  \  / 

Vl 

V2 

where  ai  and  0.2  are  the  first  two  columns  of  AT^^.  Note  that  since  is  post-multiplying 
A,  the  diagonal  elements  afai  and  020:2  not  unity  in  this  case. 

Equation  (3.27)  is  the  same  form  as  the  stationary  platform  wheel  momenta  ellipsoid 
Eq.  (3.4)  for  the  case  of  two  wheels,  expressed  in  Tv  Thus,  the  same  control  law  as 
Eq.  (3.20)  applies  to  the  three  wheel  case,  with  z/3  =  0 


/  T 

a{a.2 


U  =  € 


a.'^<X2 


0 


\ 


— Q'^cKi  —afa2  0 
0  0  0 


1/  =  eBi/ 


(3.28) 


Now,  by  solving  Eq.  (3.25)  for  i/  and  recognizing  that  we  can  substitute 

into  Eq.  (3.28)  to  obtain 


(3.29) 


or 

(3.30) 

which  is  the  sub-optimal  wheel  torque  control  law  for  a  vehicle  using  three  momentum 
wheels.  Use  of  this  control  law  with  ‘‘small”  c  should  result  in  a  reorientation  from  an 
initial  to  a  final  fij  with  little  induced  platform  motion. 

The  time  required  to  accomplish  a  stationary  platform  maneuver  is  certainly  of  in¬ 
terest  to  us,  as  this  will  be  a  significant  factor  in  assessing  the  utility  of  these  control  laws 


/i  =  eT^vBTl^fi 
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for  real  applications.  As  with  the  direct  trajectory,  we  assume  to  =  0  and  that  we  are  given 
initial  and  final  wheel  momenta  and  fi ^ . 

We  return  to  the  wheel  torque  control  law  in  the  frame,  given  in  Eq.  (3.28), 
eliminating  the  state  1/3  since  it  is  zero  and  constant.  This  results  in  a  system  of  two  linear 
constant  co-efficient  equations.  The  solution  (at  time  tj)  to  this  linear  equation  is 

vj  =  e^^^^Uo  =  (3.31) 


where  the  columns  of  V  are  the  eigenvectors  of  B,  is  a  diagonal  matrix  with  diagonal 

elements  and  A,  are  the  eigenvalues  of  B.  Putting  B  into  the  form 


B  = 


Bn  Bi2 

-B21  —Bn 


(3.32) 


allows  us  to  see  that  the  eigenvalues  will  be  purely  imaginary  as  long  as  B12B21  >  B^^. 
In  the  case  where  the  off-diagonal  elements  are  1  and  the  diagonal  elements 

are  always  less  than  1  as  long  as  the  wheels  are  non-colinear.  In  the  more  general  case, 
the  elements  Bij  are  formed  from  inner  products  of  the  first  two  columns  of  AT^^,  where 
is  orthonormal,  and  we  find  that  in  all  the  cases  presented  herein  the  eigenvalues  are 
indeed  imaginary. 

Since  the  eigenvalues  of  B  are  purely  imaginary  (±j B12B21  —  B^^),  we  note  that 
is  a  rotation  matrix  [21:212].  If  we  define  0  =  €{y/Bi2B2i  -  Bl^)tf,  and  recognize 
that  by  definition  of  that  i/„  =  (1  0)^,  then  Eq.  (3.31)  can  be  written 


cos  (f)  sin  (f) 
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Denoting  the  elements  of  the  2x2  matrix  V  ^  by  V,j,  and  the  matrix  resulting  from  V 
by  V  =  (ui  ^2)^,  Eq.  (3.33)  can  be  reduced  through  some  manipulation  to  the  form 
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which  can  be  solved  for  (j)  to  obtain 


^21^1  +  Vll^2 

Vii^l  +  V2\V2 
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(3.35) 


Recalling  the  definition  of  (j)  above  and  the  relationship  between  t  and  t  from  Eq.(2.27), 
this  can  be  solved  for  the  dimensional  final  time 
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The  magnitude  of  the  small  torque  parameter  6  vs.  if  is  obviously  found  directly  from  this 
equation,  and  is  provided  here  to  motivate  discussion. 
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Using  this  approach,  to  obtain  positive  final  time  values,  the  small  parameter  e  must 
be  negative. 

Fig.  3.16  is  a  plot  of  |6|  vs.  tj  for  the  Hubble  Space  Telescope  using  the  sub-optimal 
control  law.  The  boundary  conditions  and  vehicle  parameters  are  the  same  as  for  the  direct 
maneuver  depicted  in  Fig.  3.14,  with  4  ~  1-8  X  10^  kg  m?,  using  wheels  >Vi  -  FVa  with 
=  [1  0  0],  /Xy  =  [0  1  0]  for  both  ho  =  100  and  200  kg  rri^/s.  This  corresponds  to  a  90® 
maneuver,  so  that  for  e  =  0.1,  the  spacecraft  rotates  at  approximately  0.009®/^. 

As  with  the  direct  trajectory,  we  note  a  significant  trade-off  between  final  time  and 
torque  magnitude,  and  again  see  that  reducing  Id  ho  and/or  increasing  tj  reduces  the 
required  torques  for  this  maneuver. 

It  might  be  tempting  to  conclude  that  since  the  direct  trajectories  have  shorter  path 
lengths  in  than  sub-optimal  trajectories  for  identical  initial  and  final  conditions,  the 
maneuver  times  would  also  be  less.  However,  as  Figs.  3.17  and  3.18  show,  this  is  not  always 
the  case.  Both  figures  are  for  the  Hubble  Space  Telescope  using  wheels  Wi  —Ws,  with  equal 
total  angular  momentum  and  initial  conditions  (/x^  =  [10  0]^).  However,  whereas  the  final 
condition  in  the  first  case  in  /x^  =  [0  1  0]^,  in  the  second  it  is  /x^  ^  [0.715  —1.086  —0.074]^. 
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Figure  3.16  Sub-optimal  Trajectory  -  Small  Torque  Parameter  e  vs,  ij  for  varying  ho 

The  explanation  for  this  lies  in  the  fact  that  the  average  wheel  momenta  for  sub-optimal 
trajectories  are  generally  larger  than  for  direct  trajectories,  and  can  thus  be  traversed 
more  quickly.  As  demonstrated,  however,  the  relationship  between  e,  f/,  and  the  control 
trajectory  depends  on  the  momentum  wheels  used  and  the  boundary  conditions. 

Now  that  we  have  developed  the  wheel  torque  control  laws  for  both  direct  and  sub- 
optimal  trajectories,  let  us  run  some  example  simulations  for  both  the  GPS  Block  HR  and 
Hubble  Space  Telescope  vehicles,  in  an  environment  free  from  external  perturbing  torques. 
This  will  allow  us  to  compare  both  types  of  maneuvers,  and  particularly  to  assess  the 
utility  of  the  sub-optimal  control  law  for  “real”  vehicles.  Then,  in  Chapter  5,  external 
perturbing  torques  will  be  introduced  to  allow  assessment  of  the  sub-optimal  control  law 
for  “real”  vehicles  in  a  more  “realistic”  environment. 
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IV.  Unperturbed  Maneuvers 


In  this  chapter,  specific  examples  of  maneuvers  using  the  GPS  Block  HR  and  Hubble 
Space  Telescope  vehicles  will  be  presented.  These  two  vehicles  were  chosen  because  they 
both  employ  at  least  three  momentum  wheels  for  attitude  control,  but  in  varying  degrees. 
While  momentum  wheels  are  used  only  to  offset  small  perturbative  torques  on  the  GPS 
Block  HR  vehicle,  they  are  the  primary  attitude  control  mechanism  for  the  Hubble  Space 
Telescope.  The  momentum  storage  capacities  of  the  vehicles’  momentum  wheels  are  thus 
quite  different,  and  the  vehicles  can  be  seen  to  represent  two  ends  of  the  spectrum  for 
momentum  storage  device  use. 

Appendix  D  provides  a  functional  block  diagram  of  the  simulation  computer  program 
used  to  generate  the  results  in  this  section.  The  main  program  (rotor_m.m)  first  calls  the 
data  definition  file  (dat-gps.m,  datJiub.m)  to  load  vehicle  physical  parameters,  orbital 
parameters,  and  simulation  control  parameters.  Next,  after  calculating  the  appropriate 
wheel  torque  control  law  for  the  maneuver  (direct  or  sub-optimal),  the  program  sends  the 
appropriate  equations  of  motion  (with  or  without  perturbing  forces)  to  the  integration 
subroutine  (ode45.m).  After  integration  is  complete,  important  data  is  output  in  either 
text  or  graphical  format. 

The  simulation  program  has  the  capability  to  output  more  than  10  different  data 
plots,  in  both  dimensional  and  non-dimensional  values.  However,  the  most  important  plots 
that  will  be  used  in  our  analysis  are  the  following;  Torque  trajectory  relative  to  stationary 
platform  wheel  momenta  ellipsoid,  torque  trajectory  relative  to  spacecraft  total  momentum 
sphere,  momentum  wheel  angular  velocities  vs.  time,  spacecraft  angular  velocities  vs.  time, 
and  three  dimensional  spacecraft  attitude  plots. 

To  minimize  the  number  of  examples  that  will  be  presented,  let  us  first  state  the 
objectives  we  wish  to  meet  in  the  following  two  sections.  For  GPS,  we  will  first  perform 
a  maneuver  from  a  specific  stationary  attitude  to  another  (or  close  to  it,  depending  upon 
final  state  errors),  using  both  the  direct  control  law  and  sub-optimal  control  law,  with 
identical  e  and  h„,  to  show  the  difference  in  final  state  error  and  intermediate  angular 
velocities  for  each  approach.  Next,  we  will  calculate  the  maximum  allowable  ho  for  each 
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vehicle,  then  perform  the  same  maneuver  using  the  sub-optimal  control  law  several  times, 
varying  the  magnitude  of  e,  to  identify  what  value  of  c  is  reasonable  in  terms  of  final  state 
errors,  maneuver  time,  and  intermediate  angular  velocities. 

4.1  GPS  Block  HR  Examples 

As  mentioned  previously,  the  GPS  Block  HR  vehicle  uses  expendable  fuel  thrusters 
for  large  angle  maneuvers,  and  employs  four  reaction  wheels,  arranged  in  a  pyramid  fash¬ 
ion,  to  offset  small  perturbative  torques  acting  on  the  vehicle.  The  wheels  are  augmented 
by  a  magnetic  torquing  system  which  is  used  for  momentum  dumping.  The  vehicle  physi¬ 
cal  parameters  are  detailed  in  Appendix  B.  The  difference  in  application  between  reaction 
wheels  and  momentum  wheels  (nominal  angular  rates,  commandable  torques)  will  be  ig¬ 
nored  in  the  following  simulations,  to  better  allow  comparison  between  GPS  and  Hubble 
Space  Telescope. 

We  first  note  that  the  maximum  momentum  a  single  wheel  can  generate  is  hmax  ^ 
±8  kg  m?  js.  Thus,  to  allow  arbitrary  reorientations,  ho  should  be  set  as  less  than  8  kg  m? I s 
in  the  simulations.  The  ratio  of  wheel  inertia  to  maximum  principal  inertia  of  the  vehicle 
is  ijimax  ^  0.0046.  Although  not  derived  herein,  it  can  be  shown  that  if  we  assume 
a  single  wheel  is  attempting  to  absorb  the  angular  momentum  for  a  pure  (stable)  spin 
about  the  major  axis  of  the  vehicle,  and  that  the  wheel  is  parallel  to  the  major  axis,  the 
maximum  initial  rate  of  the  vehicle  such  that  the  wheel  can  achieve  a  stationary  platform 
is  approximately  4.014  x  10“®  radfs.  Similar  analysis  for  a  minor  axis  spin  results  in  a 
maximum  angular  rate  of  approximately  1.030  X  10“^  radfs.  While  this  analysis  is  based 
on  several  simplifying  assumptions,  it  does  allow  us  to  see  that  only  small  vehicle  angular 
velocities  can  be  offset  using  the  momentum  wheels,  which  requires  that  final  state  errors 
be  within  or  below  this  range  to  be  controllable  using  the  momentum  wheels  alone. 

To  be  realistic,  we  must  also  be  concerned  about  the  magnitude  of  the  torques  used, 
since  in  actual  applications  we  cannot  exceed  the  maximum  nor  go  below  the  minimum 
commandable  values.  Using  Eq.  (2.26),  noting  the  values  for  minimum  and  maximum 
commandable  torques  from  Appendix  B,  the  known  value  of  Ic  =  4.5266  X  10^  kg  m?,  and 
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ho  =  ^  kg  iritis  results  in 


_  Smax^c  ^  Of; 

^max  —  ~  Iz.OD 


hl_ 

dmin^c 


0.05 


(4.1) 

(4.2) 


The  purpose  of  the  following  simulations  is  both  to  compare  direct  and  sub-optimal 
control  laws  for  attitude  maneuvers,  as  well  as  demonstrate  the  growth  in  final  state  errors 
with  increasing  torque  scaling  parameter  e. 

Discussion  and  comparison  of  the  simulation  results  for  both  GPS  Block  HR  and  the 
Hubble  Space  Telescope  is  deferred  to  section  4.3. 
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GPS  Example  1  -  Direct  vs.  Sub-optimal  Trajectories 

As  our  first  example,  we  will  compare  direct  and  sub-optimal  trajectory  maneuver 
times,  angular  velocities,  and  final  state  errors  for  a  combined  45°  rotation  about  the  ba  and 
bi  axis  with  Wi  initially  containing  all  the  angular  momentum,  using  wheels  Wi,  W2, 
Table  4.1  summarizes  the  simulation  parameters.  The  results  are  summarized  in  Figs. 
4.1-4.12  and  Table  4.2. 

Table  4.1  GPS  Example  1  -  Direct  vs.  Sub-optimal  Simulation  Parameters 


parameter 

value 

A 

(ai  a2  a4) 

h/Q 

5  kg  m?  Js 

e 

0.01 

0,0,0 

45°,45°,0 

(1  0  0) 

Table  4.2  GPS  Example  1  -  Direct  vs.  Sub-optimal  Simulation  Results 


result 

Direct 

Sub-optimal 

h 

|a;| 

1  \max 

\Cb»  1 

1  ^imax 

24  hr  22  min 

3.461  X  10"^  rad/s 
590.3  rad/s 

2.788  X  10“®  rad/s 
-53.30,119.48,-5.133° 

28  hr  38  min 
1.895  X  10"®  rad/s 
613.2  rad/s 

1.526  X  10-®  rad/s 
48.35,44.04,-2.27° 
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Figure  4.1  GPS  Example  1  -  (Direct) 
Torque  Trajectory  on  Wheel 
Momenta  Ellipsoid 
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Figure  4.3  GPS  Example  1  -  (Direct) 
Torque  Trajectory  on  Total 
Momentum  Sphere 


Figure  4.2  GPS  Example  1  -  (Sub-opt) 
Torque  Trajectory  on  Wheel 
Momenta  Ellipsoid 
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Figure  4.4  GPS  Example  1  -  (Suh-opt) 
Torque  Trajectory  on  Total 
Momentum  Sphere 
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Figure  4.5  GPS  Example  1  -  (Direct) 
Wheel  Angular  Velocities 


Figure  4.7  GPS  Example  1  -  (Direct) 
Vehicle  Angular  Velocities 
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Figure  4.6  GPS  Example  1  -  (Sub-opt) 
Wheel  Angular  Velocities 


Figure  4.8  GPS  Example  1  -  (Sub-opt) 
Vehicle  Angular  Velocities 
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Figure  4.9  GPS  Example  1  -  Initial 
Attitude 


Figure  4.10  GPS  Example  1  -  Desired 
Final  Attitude 


Figure  4.11  GPS  Example  1  -  (Direct) 
Final  Attitude 


Figure  4.12  GPS  Example  1  -  (Sub^opt) 
Final  Attitude 


GPS  Example  2  -  Sub-optimal  Trajectory 

In  the  second  example,  we  will  now  address  the  sub-optimal  trajectory  only,  this  time 
increasing  e  to  0.05.  The  other  simulation  parameters  are  identical  to  the  first  example. 
Table  4.3  summarizes  the  simulation  parameters.  The  results  are  summarized  in  Figs. 
4.13-4.16  and  Table  4.4. 

Table  4.3  GPS  Example  2  -  Sub-optimal  Simulation 
Parameters 


parameter 

value 

A 

Hq 

€ 

4^0  ?  ^0?  "00 

Mo 

(aj  a.2  a4) 

5  kg  m^/s 
0.05 

0,0,0 
45“,45®,0 
(1  0  0) 

Table  4.4  GPS  Example  2  -  Sub-optimal  Simulation  Results 


result 

Sub-optimal 

h 

|d;| 

1  \max 

Imax 

5  hr  43  min 
9.482  X  10-®  rad/s 
613.2  rad/s 

7.172  X  10"®  radfs 
49.12,43.96,-3.84“ 
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ha  (kg*m"2/s) 


Figure  4,13  GPS  Example  2  -  (Sub¬ 
opt)  Torque  Trajectory  on 
Wheel  Momenta  Ellipsoid 


Figure  4.14  GPS  Example  2  -  (Sub-opt) 
Torque  Trajectory  on  Total 
Momentum  Sphere 


Figure  4.15  GPS  Example  2  -  (Sub-opt) 
Wheel  Angular  Velocities 


Figure  4.16  GPS  Example  2  -  (Sub-opt) 
Vehicle  Angular  Velocities 


GPS  Example  3  -  Sub-optimal  Trajectory 


In  the  last  example  for  GPS  Blk  HR,  we  again  address  the  sub-optimal  trajectory 
only,  this  time  increasing  e  to  0.10.  The  other  simulation  parameters  are  identical  to  the 
first  two  examples.  Table  4.5  summarizes  the  simulation  parameters.  The  results  are 
summarized  in  Figs.  4.17-4.20  ajid  Table  4.6. 


Table  4.5  GPS  Example  3  -  Sub-optimal  Simulation 
Parameters 


parameter 

value 

A 

Hq 

6 

4^0  y  ^0?  4^0 

(ai  b.2  a4) 

5  kg  m? /s 
0.10 

0,0,0 
45°,45“,0 
(1  0  0) 

Table  4.6  GPS  Example  3  -  Sub-optimal  Simulation  Results 


result 

Sub-optimal 

h 

\Cb\ 

\  \max 

1 

1  ^\max 

2  hr  52  min 
1.897  X  10-4 

613.2  rad/s 

7.999  X  10“®  rad/s 
48.58,43.34,-3.64° 
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ws  (rad^sec) 


Figure  4.17  GPS  Example  3  -  (Sub¬ 
opt)  Torque  Trajectory  on 
Wheel  Momenta  Ellipsoid 


Time  (seconds) 

Figure  4,19  GPS  Example  3  -  (Sub-opt) 
Wheel  Angular  Velocities 


Figure  4.18  GPS  Example  3  -  (Sub-opt) 
Torque  Trajectory  on  Total 
Momentum  Sphere 


Time  (seconds) 

Figure  4.20  GPS  Example  3  -  (Sub-opt) 
Vehicle  Angular  Velocities 


J^.2  Hubble  Space  Telescope  Examples 

Unlike  GPS,  the  Hubble  Space  Telescope  uses  momentum  wheels  exclusively  for 
attitude  control,  both  for  large  angle  maneuvers  and  fine  pointing  control.  Although  data 
on  the  wheel  inertias  and  angular  rates  were  available,  their  orientation  within  the  vehicle 
was  not.  Thus,  the  same  arrangement  as  for  GPS  Block  IIR  was  assumed  for  simplicity 
(see  Appendix  C).  A  magnetic  torquing  system  is  used  for  momentum  dumping  [2]. 

The  Hubble  pointing  control  system  is  currently  able  to  slew  the  vehicle  through  a 
90°  reorientation  using  and  eigen-axis  maneuvering  scheme  in  approximately  18  minutes 
[2],  with  0.01  arc-second  pointing  accuracy.  It  will  thus  be  of  interest  to  compare  the 
sub-optimal  control  results  to  these  values. 

The  maximum  momentum  a  single  wheel  can  generate  is  hmax  ^  ±263.4  kg  m? j s. 
While  this  is  obviously  much  greater  on  an  absolute  scale,  the  ratio  of  wheel  inertia  to 
maximum  principal  inertia  of  the  vehicle  is  actually  smaller  than  for  GPS  Block  IIR  at 
islimax  ~  0.0033.  To  allow  arbitrary  reorientations,  hg  should  be  set  somewhat  less  than 
263.4  kg  m? Js  in  the  simulations.  In  all  the  examples  below,  hg  is  set  at  200  kg  rn^/s. 
As  we  shall  see  in  Example  3,  due  to  the  wheels  chosen,  this  results  in  there  being  small 
regions  of  the  wheel  momenta  ellipsoid  that  are  unachievable,  resulting  in  excessively  high 
wheel  speeds  for  trajectories  that  traverse  these  regions.  The  unachievable  regions  are  only 
shown  in  the  figures  for  Example  3,  but  they  are  present  for  Examples  1  and  2  as  well. 

Although  not  derived  explicitly,  it  can  be  shown  that  if  we  assume  a  single  wheel 
is  attempting  to  absorb  the  angular  momentum  for  a  pure  (stable)  spin  about  the  major 
axis  of  the  vehicle,  and  that  the  wheel  is  parallel  to  the  major  axis,  the  maximum  initial 
rate  of  the  vehicle  such  that  the  wheel  can  achieve  a  stationary  platform  is  approximately 
3.299  X 10“^  radfs.  Similar  analysis  for  a  minor  axis  spin  results  in  a  maximum  angular  rate 
of  approximately  1.215  X  10“^  radjs.  While  this  analysis  is  based  on  several  simplifying 
assumptions,  it  does  allow  us  to  see  that  only  small  vehicle  angular  velocities  can  be  offset 
using  the  momentum  wheels,  which  requires  that  final  state  errors  be  within  or  below  this 
range  to  be  controllable  using  the  momentum  wheels  alone. 
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To  be  realistic,  we  must  also  be  concerned  about  the  magnitude  of  the  torques  used, 
since  in  actual  applications  we  cannot  exceed  the  maximum  nor  go  below  the  minimum 
commandable  values.  Using  Eq.  (2.26),  noting  the  assumed  values  for  minimum  and  maxi¬ 
mum  commandable  torques  from  Appendix  C,  the  known  value  of  4  =  1.7815  X  10®  kg  m^, 
and  ho  =  200  kg  m?/s  results  in 


^max 


^min 


«  0.891 

(4.3) 

hl_ 

«  0.001 

(4.4) 

The  purpose  of  the  following  simulations  is  both  to  compare  direct  and  sub-optimal 
control  laws  for  attitude  maneuvers,  as  well  as  demonstrate  the  growth  in  final  state  errors 
with  increasing  torque  scaling  parameter  e.  The  last  example  is  also  useful  in  comparing  the 
sub-optimal  control  performance  to  the  eigen-axis  maneuvers  mentioned  above.  Discussion 
and  comparison  of  the  simulation  results  for  both  GPS  Block  HR  and  the  Hubble  Space 
Telescope  is  deferred  to  section  4.3. 
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Hubble  Example  1  -  Sub-optimal  Trajectory 


In  the  first  example  for  Hubble,  we  will  address  the  sub-optimal  trajectory  only  for 
a  combined  45°  rotation  about  the  ba  and  bi  axes,  with  e  set  at  0.01.  Momentum  wheels 
Wi  -  Ws  are  used,  with  Wi  initially  containing  all  the  momentum  for  the  vehicle.  Table 
4.7  summarizes  the  simulation  parameters.  The  results  are  summarized  in  Figs.  4.21-4.26 
and  Table  4.8. 


Table  4.7  Hubble  Example  1  -  Sub-optimal  Simulation 
Parameters 


parameter 

value 

A 

hjQ 

c 

Mo 

(ai  as) 

200  kg  m?  js 
0.01 

0,0,0 
45°,45°,0 
(1  0  0) 

Table  4.8  Hubble  Example  1  -  Sub-optimal  Simulation 
Results 


result 

Sub-optimal 

h 

\Cj\ 

1  \max 

IcD*  1 

1  ^  \max 

Hh)Lar 

45  hr  5  min 

1.506  X  10"®  rad/s 
287.4  radjs 

8.871  X  10"®  radjs 
72.89,37.74,-18.44° 
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ws  (rad/sec) 


Figure  4.21  Hubble  Example  1  -  (Sub¬ 
opt)  Torque  Trajectory  on 
Wheel  Momenta  Ellipsoid 


Time  (seconds) 


Figure  4.23  Hubble  Example  1  -  (Sub¬ 
opt)  Wheel  Angular 
Velocities 


Figure  4.22  Hubble  Example  1  -  (Sub¬ 
opt)  Torque  Trajectory  on 
Total  Momentum  Sphere 


Time  (seconds) 

Figure  4.24  Hubble  Example  1  -  (Sub¬ 
opt)  Vehicle  Angular 
Velocities 


Hubble  Example  2  -  Sub-optimal  Trajectory 


In  the  next  example,  we  again  address  the  sub-optimal  trajectory  only,  this  time 
increasing  e  to  0.5.  The  other  simulation  parameters  are  identical  to  the  first  example. 
Table  4.9  summarizes  the  simulation  parameters.  The  results  are  summarized  in  Figs. 
4.27-4.32  and  Table  4.10. 


Table  4.9  Hubble  Example  2  -  Sub-optimal  Simulation 
Parameters 


parameter 

value 

A 

ho 

€ 

•^0)  ^0) 

Mo 

(ai  a.2  as) 
200  kg  m?js 
0.5 

0,0,0 
45%45°,0 
(1  0  0) 

Table  4.10  Hubble  Example  2  -  Sub-optimal  Simulation 
Results 


result 

Sub-optimal 

h 

Ifimax 

Imar 

Hh)La. 

54  min 

6.909  X  10"^  radjs 
287.4  radjs 

3.740  X  10"^  radjs 
72.94,38.25,-2.89° 
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ws  (rad/sec) 


Figure  4.27  Hubble  Example  2  -  (Sub¬ 
opt)  Torque  Trajectory  on 
Wheel  Momenta  Ellipsoid 


Figure  4.29  Hubble  Example  2  -  (Sub¬ 
opt)  Wheel  Angular 
Velocities 


Figure  4.28  Hubble  Example  2  -  (Sub¬ 
opt)  Torque  Trajectory  on 
Total  Momentum  Sphere 


xIO"* 


Figure  4.30  Hubble  Example  2  -  (Sub¬ 
opt)  Vehicle  Angular 
Velocities 


le  (sec) 


9 


Hubble  Example  3  -  Direct  vs.  Sub-optimal  Trajectories 


In  this  example  for  Hubble,  both  the  direct  and  the  sub-optimal  control  laws  are  used 
to  attempt  a  positive  90°  “slewing”  maneuver  about  the  bi  axis.  To  keep  both  maneuver 
times  and  state  errors  reasonable,  e  is  set  to  0.1.  Wheels  >Vi  -  Ws  are  again  used  to 
perform  the  maneuver,  with  W2  initially  containing  all  the  angular  momentum  for  the 
vehicle.  Table  4.11  summarizes  the  simulation  parameters.  The  results  are  summarized  in 
Figs.  4.33-4.44  and  Table  4.12. 


Table  4.11  Hubble  example  3  -  Direct  vs.  Sub-optimal  Sim¬ 
ulation  Parameters 


parameter 

value 

A 

(ai  a2  aa) 

h/Q 

200  kg  m^/s 

e 

0.10 

^01  ^07  '4^0 

0,0,0 

0,90",0 

(0  1  0) 

Table  4.12  Hubble  example  3  -  Direct  vs.  Sub-optimal  Sim¬ 
ulation  Results 


result 

Sub-optimal 

Direct 

h 

IcDl 

1  \max 

1^^  Imar 

Hh)La. 

5  hr  30  min 
1.702  X  10-4 

336.7  rad/s 
7.042  X  10"®  rad/s 
-0.65,89.39,0.28” 

3  hr  30  min 

7.625  X  10-4 

238.1  rad/s 

2.050  X  10-4 
-22.69,96.15,16.49° 
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Figure  4.33  Hubble  Example  3  -  (Sub¬ 
opt)  Torque  Trajectory  on 
Wheel  Momenta  Ellipsoid 

x3 


Figure  4.34  Hubble  Example  3  -  (Di¬ 
rect)  Torque  Trajectory  on 
Wheel  Momenta  Ellipsoid 

X3 


Figure  4.35  Hubble  Example  3  -  (Sub-  Figure  4.36  Hubble  Example  3  -  (Di- 

opt)  Torque  Trajectory  on  rect)  Torque  Trajectory  on 

Total  Momentum  Sphere  Total  Momentum  Sphere 
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Figure  4.37  Hubble  Example  3  -  (Sub¬ 
opt)  Wheel  Angular 
Velocities 


Figure  4.39  Hubble  Example  3  -  (Sub¬ 
opt)  Vehicle  Angular 
Velocities 


Figure  4.38  Hubble  Example 

3  -  (Direct)  Wheel  Angular 
Velocities 


Figure  4.40  Hubble  Example  3 
-  (Direct)  Vehicle  Angular 
Velocities 
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Figure  4.41  Hubble  Example  S  -  Initial  Figure  4.42  Hubble  Example  3  -  Desired 

Attitude  Final  Attitude 


Figure  4.43  Hubble  Example  3  -  (Sub-  Figure  4.44  Hubble  Example  3  -  (Di- 

opt)  Final  Attitude  rect)  Final  Attitude 


Hubble  Example  4  -  Sub-optimal  Trajectory 


In  the  final  example  for  Hubble,  we  will  perform  the  same  maneuver  as  in  Example  3, 
this  time  using  wheels  >Vi ,  W2,  and  W4,  with  W2  initially  containing  all  the  momentum  for 
the  vehicle.  Table  4.13  summarizes  the  simulation  parameters.  The  results  are  summarized 
in  Figs.  4.45-4.50  and  Table  4.14. 


Table  4.13  Hubble  Example  4  -  Sub-optimal  Simulation 
Parameters 


parameter 

value 

A 

Hq 

€ 

4^0  9  ^09  4^0 

</>/9^/9^/ 

Mo 

(ai  a2  34) 
200  kg  rn^js 
0.1 

0,0,0 
0,90°,0 
(0  1  0) 

Table  4.14  Hubble  Example  4  -  Sub-optimal  Simulation 
Results 


result 

Sub-optimal 

h 

|a;| 

1  imax 

1^*  1 

1  ^Imaar 

\  \  J  /  \max 

3  hr  53  min 
2.406  X  10"^  rad/s 
238.1  rad/s 
4.816  X  10"®  rad/s 
-0.82,90.43,0.70“ 
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mu3 


X3 


Figure  4.45  Hubble  Example  4  -  (Sub¬ 
opt)  Torque  Trajectory  on 
Wheel  Momenta  Ellipsoid 


Figure  4.46  Hubble  Example  4  "  (Sub¬ 
opt)  Torque  Trajectory  on 
Total  Momentum  Sphere 


Time  (seconds) 


xIO*^ 


Figure  4.47  Hubble  Example  4  -  (Sub¬ 
opt)  Wheel  Angular 
Velocities 


Figure  4.48  Hubble  Example  4  -  (Sub¬ 
opt)  Vehicle  Angular 
Velocities 
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Figure  4.50  Hubble  Example  4  -  (Sub-opt)  Final  Attitude 
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4-3  Summary  of  Results 

GPS  Example  1  and  Hubble  Example  3  both  compared  direct  and  sub-optimal  ma¬ 
neuvers  for  a  specified  initial  conditions  and  desired  final  wheel  momenta  and  vehicle  atti¬ 
tude  (See  Tables  4.2  and  4.12).  As  discussed  in  Section  3.5,  the  relative  time  to  complete  a 
maneuver  using  both  approaches  varies  depending  on  the  parameters  of  the  problem.  For 
these  two  somewhat  dissimilar  examples,  the  sub-optimal  trajectory  required  between  2 
and  4  hours  longer  to  complete  the  maneuver  than  the  direct  trajectory. 

For  both  examples,  the  the  maximum  wheel  speeds  for  the  direct  trajectory  were  well 
below  the  maximum  allowable  values,  with  the  sub-optimal  trajectory  requiring  slightly 
higher  speeds.  While  the  wheel  speeds  for  GPS  were  well  below  the  943  radfs  limit  for  both 
approaches,  the  sub-optimal  trajectory  for  the  Hubble  example  exceeded  the  maximum  of 
314  rad/s  for  a  small  portion  of  the  simulation.  A  look  at  the  wheel  speed  histories  in 
Figures  4.5, 4.6,  4.37,  and  4.38  shows  that  it  is  due  to  the  quadratic  form  of  the  sub-optimal 
control  law  that  the  maximum  is  reached  in  the  sub-optimal  cases.  Thus,  the  sub-optimal 
control  law  is  infeasible  due  to  wheel  speed  limitations  for  smaller  values  of  total  angular 
momentum  than  is  the  direct  trajectory.  Recalling  our  discussion  of  achievable  momenta 
ellipsoid  regions  in  Section  3.2,  it  is  apparent  that  trajectories  that  have  initial  and  final 
states  near  the  axes  planes  in  and  near  the  “ends”  of  the  major  axis  of  the  ellipsoid  will 
traverse  regions  of  the  ellipsoid  requiring  larger  wheel  momenta  due  to  the  nature  of  the 
sub-optimal  control  law.  In  such  cases,  other  trajectories  on  the  ellipsoid  surface  would  be 
necessary,  or  different  sets  of  momentum  wheels  could  be  used.  Hubble  Example  4  uses  >Vi, 
W2,  and  yVi  to  perform  the  same  maneuver  with  the  same  total  angular  momentum,  using 
the  sub-optimal  control  law.  Note  that  using  these  three  wheels,  the  sub-optimal  trajectory 
no  longer  traverses  the  unachievable  regions  of  the  ellipsoid,  the  final  attitude  is  very  close 
to  desired,  and  the  resulting  maximum  wheel  speeds  are  below  the  limit  throughout  the 
simulation  (Table  4.14).  This  clearly  demonstrates  the  advantages  of  having  a  redundant 
non-coplanar  wheel  in  momentum  limited  situations  such  as  Example  3. 

In  terms  of  maximum  vehicle  angular  velocities  during  maneuvers,  the  sub-optimal 
trajectory  is  clearly  superior.  While  absolute  angular  velocities  remained  small  for  both 
direct  and  sub-optimal  trajectories  (due  to  small  wheel  torques),  the  relative  magnitude  for 
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the  direct  trajectories  was  between  4  and  18  times  larger.  This  is  certainly  not  surprising 
given  the  definition  of  the  stationary  platform  condition  and  the  way  in  which  each  torque 
trajectory  traverses  the  stationary  platform  wheel  momenta  ellipsoid. 

What  is  interesting  are  the  final  vehicle  angular  velocities  using  each  approach  for 
the  two  vehicles.  One  might  expect  that  the  sub-optimal  trajectories,  since  they  maintain 
smaller  angular  velocities  throughout  the  maneuver,  would  “end  up”  with  smaller  final 
angular  velocities.  However,  a  glance  at  the  results  shows  that  this  is  not  true  in  GPS 
Example  1.  Deviations  from  the  stationary  platform  condition  (requiring  x  =  A/t  and 
w  =  0)  can  be  visualized  by  the  magnitude  of  the  “humps”  on  the  torque  trajectory  traced 
on  the  total  momentum  sphere  (Figures  4.3,4.4  and  4.35,4.36).  Thus,  larger  deviations 
(caused  by  larger  e  and/or  straying  farther  off  the  wheel  momenta  ellipsoid  surface)  result 
in  larger  humps  and  larger  variations  in  angular  velocities.  The  final  angular  velocities 
(which  can  be  thought  of  as  state  errors  from  the  nominal  =  0)  depend  on  the  distance 
between  the  actual  final  state  Xact  vs.  the  desired  final  state  x/  which  is  calculated  from 
the  stationary  platform  condition  x/  =  Afij  in  x  space  at  the  time  when  the  wheel  torques 
are  turned  off.  This  is  difficult  to  see  in  the  GPS  example  (since  e  is  much  smaller),  but  in 
the  Hubble  example  (Fig.  4.36),  it  is  obvious  that  x^ct  is  farther  from  X/  (denoted  by  the 
X  on  the  sphere)  in  the  direct  trajectory  case.  However,  because  the  period  of  the  humps 
plays  as  important  a  role  as  the  amplitude,  there  may  be  cases  where  the  final  state  just 
happens  to  be  closer  to  the  desired  state  for  the  direct  trajectory,  even  though  absolute 
variations  in  angular  velocities  are  greater.  This  situation  occurs  in  Hubble  Example  3. 

Perhaps  the  most  important  results  to  consider  are  the  errors  in  the  final  attitude 
after  each  maneuver.  As  shown  in  Table  4.12  and  Figures  4.43  and  4.44,  the  sub-optimal 
control  law  results  in  final  attitude  within  0.65®  of  desired,  while  the  direct  trajectory 
leaves  the  vehicle  over  20°  off  in  some  directions.  As  discussed  at  the  end  of  Section  3.1, 
the  final  attitude  for  both  trajectories  is  off  by  a  finite  rotation  about  the  total  angular 
momentum  vector  x,  recognizing  that  there  are  final  state  errors  as  well.  It  is  clear  that  for 
this  example,  the  sub-optimal  control  results  in  a  much  smaller  necessary  rotation  about 
X  and  thus  significantly  smaller  final  attitude  errors. 
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Turning  now  to  a  comparison  of  the  sub-optimal  trajectories  for  varying  values  of 
e  (GPS  Examples  1-3,  Hubble  Examples  1-2),  it  is  evident  that  increasing  e  for  a  given 
ho  reduces  maneuver  times  but  results  in  larger  intermediate  and  final  angular  velocities 
and  greater  final  attitude  errors.  These  results  are  certainly  not  unexpected  given  the 
discussion  earlier  in  the  chapter. 

What  is  interesting,  though,  is  that  when  e  is  increased  significantly  (e.g.  0.5  in  Hub¬ 
ble  Example  2),  maneuver  times  become  reasonable  (less  than  an  hour),  angular  velocities 
during  the  maneuver  increase,  yet  the  motion  of  the  vehicle  is  still  much  “smoother”  than 
the  direct  maneuver  with  much  smaller  torques  (see  Hubble  Example  3).  As  a  result, 
maneuvers  can  be  accomplished  in  realistic  time-frames  and  yet  maintain  small  platform 
angular  velocities  and  achieve  a  final  attitude  near  the  desired  one.  Although  maneuver 
times  using  this  approach  were  still  greater  than  can  be  achieved  currently  with  the  Hubble 
using  eigen-axis  maneuvers,  the  nature  of  the  motion  is  “smooth”  and  thus  desirable.  Ob¬ 
viously,  the  maximum  value  of  e  found  previously  cannot  be  exceeded  for  a  given  ho,  since 
this  would  result  in  commanding  the  wheel  torque  motors  beyond  their  practical  limit. 
However,  it  is  apparent  that  even  at  torques  well  below  the  maximum  commandable  (as 
in  Hubble  Example  2),  the  sub-optimal  control  law  can  effect  large  angle  reorientations  in 
reasonable  times. 

In  all  the  examples  for  GPS  and  Hubble,  the  angular  velocities  at  the  end  of  the 
maneuvers  were  all  of  order  10“'’  or  smaller,  while  the  maximum  allowable  body  rates 
calculated  in  Sections  3.6.1  and  3.6.2  were  of  order  10“®  -  10“^.  Thus,  one  might  assume 
that  the  final  angular  velocities  were  acceptably  small,  and  could  be  offset  using  additional 
momentum  wheel  control.  However,  it  is  important  to  note  that  the  analyses  previously 
performed  to  calculate  these  maximum  body  rates  were  approximate,  and  also  assumed 
that  the  wheels  were  initially  not  spinning  and  were  torqued  to  absorb  this  momentum. 
In  all  the  examples  above,  the  wheels  already  contain  some  angular  momentum  (i.e.  they 
are  spinning),  and  thus  must  absorb  this  additional  angular  momentum  due  to  the  body 
angular  rates.  Thus,  depending  on  the  magnitude  of  the  body  rates  and  the  final  wheel 
momenta,  the  additional  residual  angular  momentum  may  or  may  not  be  absorbable  by 
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the  wheels  without  momentum  dumping  (using  magnetic  torquers  and/or  thrusters  in  the 
case  of  GPS). 

Nonetheless,  the  examples  show  that  the  sub-optimal  control  is  easily  calculable,  is 
superior  to  the  constant  torque  case  in  terms  of  body  angular  rates  and  final  state  errors, 
and  can  be  used  to  achieve  large  angle  reorientations  in  reasonable  times. 

Now,  let  us  introduce  solar  pressure  and  gravity  gradient  torques  to  increase  the 
“realism”  of  the  simulations,  and  assess  the  affect  on  intermediate  angular  velocities  and 
final  state  errors. 
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V.  Perturbed  Stationary  Platform  Maneuvers 

As  discussed  during  the  development  of  the  external  perturbing  torque  equations  in 
Chapter  2,  the  magnitudes  of  solar  pressure  and  gravity  gradient  torques  vary  depending 
on  vehicle  configuration,  orbital  altitude,  and  spacecraft  attitude.  While  gravity  gradient 
torques  decrease  in  absolute  magnitude  with  altitude,  solar  pressure  torques  remain  essen¬ 
tially  constant  (since  P,u„  =  const).  Thus,  for  a  given  vehicle,  as  orbital  altitude  increases, 
the  relative  contribution  of  solar  pressure  to  total  external  torques  increases. 

Both  of  these  vehicles  are  in  nearly  circular  orbits  and  are  asymmetric  in  both  physical 
shape  and  inertia  properties.  Recalling  the  discussion  about  gravity  gradient  stability  from 
Section  2.3.1,  we  expect  gravity  to  attempt  to  reorient  both  vehicles  such  that  their  minor 
axes  are  aligned  with  the  local  vertical.  Of  course,  solar  pressure  is  also  affecting  the 
attitude  throughout  the  orbit. 

GPS  Block  HR  vehicles  are  at  much  higher  altitudes  than  Hubble  and  are  significantly 
less  massive,  and  as  a  result,  gravity  gradient  torques  are  of  roughly  the  same  magnitude 
as  solar  pressure  torques.  Using  Eq.  (2.46)  and  (2.49)  for  GPS  Block  HR  an  altitude 
of  26,560  km,  assuming  Pt  is  the  principal  frame  (this  is  nearly  true),  with  the  bi  axis 
rotated  -f  45“  in  the  ei  -  63  plane  initially,  the  solar  pressure  and  gravity  gradient  torques 
are  both  found  to  be  on  the  order  of  10"®  N  m.  However,  for  Hubble,  the  same  analysis 
at  an  altitude  of  6, 993  km  results  in  gravity  gradient  being  of  order  10“^  JV  m  at  times, 
while  solar  pressure  torques  remain  of  order  10"®  N  m.  The  solar  pressure  torques  are 
larger  on  average  for  Hubble  due  to  its  larger  surface  area. 

As  a  result  of  this  difference  in  both  relative  and  absolute  external  torques,  we  expect 
that  Hubble  will  experience  greater  fluctuation  in  angular  momentum  and  thus  larger  state 
errors  during  stationary  platform  maneuvers  than  GPS  for  equivalent  maneuver  times. 

GPS  Block  HR  vehicles  also  have  significantly  longer  orbital  periods  than  Hubble 
(w  12  hrs  vs.  w  1.6  hrs).  If  the  Earth’s  gravitational  field  were  modeled  more  accurately 
(not  uniform),  and  if  eclipse  effects  were  included,  both  solar  pressure  and  gravity  gradient 
torques  would  be  more  pronounced  for  the  Hubble  since  multiple  revolutions  could  occur 
during  a  2  hour  maneuver. 
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All  of  these  factors  result  in  net  changes  in  total  angular  momentum  (h)  of  the 
vehicle.  The  sub-optimal  wheel  torque  control  law  developed  in  Chapter  3  assumed  that  h 
was  constant,  since  it  was  based  on  the  stationary  platform  condition  x  =  A/x  which  was 
non-dimensionalized  by  an  initial  “nominal”  angular  momentum  (ho).  Changes  in  total 
angular  momentum  from  the  “nominal”  thus  result  in  a  change  in  the  “size”  of  the  total 
momentum  sphere  and  the  stationary  platform  wheel  momenta  ellipsoid.  Accordingly,  even 
though  the  sub-optimal  control  law  is  closed-loop  in  the  sense  that  the  control  is  calculated 
based  on  vehicle  states  (/x  and  x),  it  is  based  on  an  incomplete  dynamics  model. 

In  real  applications,  then,  as  the  wheel  momenta  ellipsoid  expands  and  contracts, 
some  provision  is  necessary  to  control  the  wheel  speeds  to  “stay”  on  the  ellipsoid  surface 
if  a  stationary  platform  is  desired.  Furthermore,  if  a  particular  stationary  attitude  is 
desired,  the  wheels  must  be  controlled  to  stay  at  a  particular  point  on  the  ellipsoid  surface. 
Finally,  if  a  stationary  platform  maneuver  is  attempted,  the  vehicle  states  would  need  to 
be  accurately  estimated  to  allow  constant  modification  of  the  sub-optimal  trajectory. 

Two  examples  are  presented  below  which  demonstrate  the  relative  and  absolute  ef¬ 
fects  of  solar  pressure  and  gravity  gradient  torques  during  a  stationary  platform  maneuver. 
The  first  example  is  the  same  maneuver  as  performed  in  GPS  Example  3  in  Chapter  3, 
with  perturbing  torques  included.  The  second  example  is  the  same  as  Hubble  Example  4, 
with  torques  included.  These  two  examples  were  chosen  for  several  reasons.  The  maneu¬ 
ver  times  are  both  reasonable  and  nearly  equal  («  1  X  10^  s),  final  attitudes  were  close  to 
desired,  and  final  angular  velocities  were  small  (10““*  —  10“®  radfs). 

Analysis  of  the  results  of  both  simulations  is  presented  in  Section  5.3. 
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5 A  Example  Using  GPS  Blk  IIR 

This  example  using  GPS  Block  IIR  is  the  same  as  GPS  Example  3  in  Chapter  4,  with 
solar  pressure  and  gravity  gradient  torques  acting  on  the  vehicle.  Table  5.1  summarizes 
the  simulation  parameters.  The  results  are  summarized  in  Figs.  5. 1-5.6  and  Table  5.2. 


Table  5.1  GPS  example  4  -  Perturbed  Sub-optimal  Simula¬ 
tion  Parameters 


parameter 

value 

A 

ho 

€ 

4^0  5  ^0  5  '^0 

Po 

(ai  a.2  a4) 

5  kg  m?  js 
0.10 

0,0,0 
45“,45“,0 
(1  0  0) 

Table  5.2  GPS  example  4  -  Perturbed  Sub-optimal  Simula¬ 
tion  Results 


result 

Sub-optimal 

h  _ 

number  of  orbital  periods 
|d;| 

\  \max 

lu)*  1 

1  ^\max 

Hh)La. 

2  hr  52  min 

0.24 

2.210  X  10"^  rad/s 
613.2  rad/s 

1.577  X  10-^  rad/s 
1.04, 77.39, 14.07° 
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(oes/pej)  SM 


Figure  5.1  GPS  Example  4  -  (Sub-opt) 
Perturbed  Torque  Trajectory 
on  Total  Momentum  Sphere 


Figure  5.3  GPS  Example  4  ~  (Sub¬ 
opt)  Perturbed  Wheel  Angu¬ 
lar  Velocities 


0  2000  4000  6000  8000  10000  12000 

Time  (seconds) 

Figure  5.2  GPS  Example  4  "  (Sub¬ 
opt)  Perturbed  Total  Angu¬ 
lar  Momentum  Magnitude 


Figure  5.4  GPS  Example  4  -  (Sub-opt) 
Perturbed  Vehicle  Angular 
Velocities 


5.2  Example  Using  Hubble  Space  Telescope 

This  example  for  Hubble  Space  Telescope  is  the  same  as  Hubble  Example  4  in  Chap¬ 
ter  4,  with  solar  pressure  and  gravity  gradient  torques  acting  on  the  vehicle.  Table  5.3 
summarizes  the  simulation  parameters.  The  results  are  summarized  in  Figs.  5.7-5.12  and 
Table  5.4. 


Table  5.3  Hubble  Example  5  -  Perturbed  Sub-optimal  Sim¬ 
ulation  Parameters 


parameter 

value 

A 

hjQ 

€ 

4^0 1  ^0?  4^0 

(ai  3.2  a4) 

200  kg  m^/s 
0.10 

0,0,0 
0,90%0 
(0  1  0) 

Table  5.4  Hubble  Example  5  -  Perturbed  Sub-optimal  Sim¬ 
ulation  Results 


result 

Sub-optimal 

h 

number  of  orbital  periods 

\ij\ 

1  \max 

1^’ \max 

Hh)Lax 

3  hr  53  min 

2.40 

1.582  X  10“^  rad/s 
238.1  rad/s 

5.690  X  10-^  rad/s 
53.84,66.31,-54.49“ 
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Figure  5.7  Hubble  Example  5  -  (Sub¬ 
opt)  Perturbed  Torque  Tra¬ 
jectory  on  Total  Momentum 
Sphere 


Figure  5.8  Hubble  Example  5  -  (Sub¬ 
opt)  Perturbed  Total  Angu¬ 
lar  Momentum  Magnitude 


Figure  5.9  Hubble  Example  5  -  (Sub- 
opt)  Perturbed  Wheel  Angu¬ 
lar  Velocities 


Figure  5.10  Hubble  Example  5  -  (Sub¬ 
opt)  Perturbed  Vehicle  An¬ 
gular  Velocities 
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time  (sec) 
1.399e+04 


Figure  5.12  GPS  Example  5  -  (Sub-opt)  Perturbed  Final  Attitude 


5.3  Summary  of  Results 

As  expected,  the  maximum  angular  velocities  during  the  perturbed  maneuvers  in¬ 
creased  for  both  vehicles  compared  to  the  unperturbed  maneuvers.  Again  as  predicted, 
the  increase  for  Hubble  was  more  dramatic  (w  6.56  times)  than  for  GPS  (w  1.16  times). 
The  maximum  final  angular  velocities  likewise  increased  more  for  Hubble  than  for  GPS 
(fa  11.81  times  vs.  Ri  1.97  times),  and  the  pointing  accuracy  at  the  completion  of  the  ma¬ 
neuver  is  arguably  better  for  GPS.  The  figures  show  that  Hubble  does  indeed  experience 
both  a  greater  absolute  increase  in  total  momentum  (Figs.  5.2  and  5.8)  and  a  greater  rel¬ 
ative  increase  (Figs.  5.1  and  5.7),  due  to  its  greater  mass,  surface  area,  and  lower  altitude 
orbit. 

The  maneuver  time  for  GPS  was  approximately  2.9  hours,  compared  to  nearly  3.9 
hours  for  the  Hubble.  This  equates  to  nearly  2.4  orbital  periods  for  Hubble  vs.  only  0.24  for 
GPS.  These  relatively  long  maneuvers  (compared  to  Hubble  Example  2  which  is  less  than 
an  hour)  were  intentionally  chosen  to  show  the  relative  affects  of  the  perturbing  torques 
for  these  two  vehicles.  The  same  affects  would  be  present,  although  less  obvious,  for  more 
“realistic”  wheel  torques  (e  fs  0.5)  as  well.  A  glance  at  Figs.  5.3  and  5.9  shows  that 
while  both  experience  general  growth  in  angular  velocities,  the  gravity  gradient  and  solar 
pressure  torques  are  more  evident  in  Hubble  Example  5,  particularly  for  As  expected, 
the  gravity  gradient  has  more  time  to  act  on  Hubble,  and  is  attempting  to  align  the  minor 
axis  with  the  local  vertical. 

Even  though  the  performance  of  the  sub-optimal  control  law  appears  to  be  better  for 
for  GPS  than  for  Hubble  in  the  presence  of  perturbing  torques,  in  terms  of  final  attitude 
errors  it  is  clearly  inadequate  in  its  present  implementation.  It  is  obvious  that  a  closed 
loop  controller  needs  to  be  developed  which  can  use  estimates  of  the  vehicle  states  (wheel 
angular  velocities,  vehicle  angular  velocities,  vehicle  attitude)  to  control  wheel  torques 
to  traverse  the  stationary  platform  torque  trajectory  even  when  the  momenta  ellipsoid 
is  expanding  and  the  angular  velocity  vector  is  moving  with  respect  to  inertial  space. 
However,  the  results  of  these  simulations  using  it  in  open-loop  form  indicate  that  it  is 
more  effective  for  vehicles  experiencing  smaller  perturbing  torques.  This  is  particularly 
true  since  increases  in  total  angular  momentum  will  result  in  the  need  for  increased  wheel 
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speeds  for  almost  all  maneuvers,  which  could  push  the  wheels  to  their  commandable  limits. 
For  the  Hubble,  such  limits  would  be  reached  more  quickly,  resulting  in  more  unachievable 
stationary  attitudes,  or  to  avoid  this,  the  need  for  more  frequent  momentum  dumping 
using  the  magnetic  torquers. 
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VI.  Conclusions  and  Recommendations 


6. 1  Conclusions 

This  report  focused  on  reorientations  of  a  rigid  spacecraft  with  three  cixisymmet- 
ric  momentum  wheels.  After  deriving  the  equations  of  motion  and  expressions  for  solar 
pressure  and  gravity  gradient  torques  in  Chapter  2,  two  conditions  involving  the  spacecraft 
angular  momentum  and  wheel  angular  momenta  were  presented  in  Chapter  3  such  that  the 
vehicle  could  be  held  stationary  with  respect  to  inertial  space  using  the  momentum  wheels 
alone.  After  deriving  a  constant  wheel  torque  control  law  (direct  trajectory)  to  transfer 
wheel  angular  momentum  without  regard  to  intermediate  angular  velocities,  a  sub-optimal 
wheel  torque  control  law  was  presented  that  kept  angular  velocities  small  throughout  the 
maneuver. 

Simulation  runs  using  both  the  direct  and  sub-optimal  control  law  were  performed 
in  Chapter  4  for  the  GPS  Block  HR  and  Hubble  Space  Telescope  vehicles,  in  an  environ¬ 
ment  free  from  external  perturbing  torques.  As  expected,  intermediate  angular  velocities, 
final  angular  velocities,  and  final  attitude  errors  were  significantly  smaller  using  the  sub- 
optimal  control  law  vs.  using  constant  torques.  Also,  as  the  magnitude  of  the  control 
torques  increased  using  the  sub-optimal  approach,  angular  velocities  (intermediate  and 
final)  increased,  as  did  the  final  attitude  error.  However,  by  increasing  the  small  torque 
parameter  e  to  0.5,  a  45°  reorientation  of  the  Hubble  Space  Telescope  about  two  axes 
was  accomplished  in  54  minutes  using  the  sub-optimal  control  law,  with  acceptably  small 
intermediate  and  final  angular  velocities  and  final  attitude  close  to  desired.  The  resulting 
motion  was  quite  “smooth”  compared  to  the  direct  maneuver,  and  showed  the  applicability 
of  this  type  of  control  scheme  in  real  applications. 

Gravity  gradient  and  solar  pressure  torques  were  introduced  in  Chapter  5,  and  two 
previous  simulations  were  re-run  with  these  torques  present  (one  for  GPS,  one  for  Hubble). 
As  the  preliminary  analysis  predicted,  the  effects  on  Hubble  were  much  greater  than  on 
GPS  Block  HR  for  equivalent  time  spans,  due  to  Hubble’s  larger  surface  area,  greater 
mass,  and  lower  altitude  orbit.  Intermediate  and  final  angular  velocities  increased  6-11 
times  over  the  unperturbed  maneuver,  compared  to  only  1-2  times  for  GPS.  Although 
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difficult  to  judge,  the  final  attitude  for  GPS  appeared  to  be  slightly  closer  to  desired  than 
Hubble. 

All  things  considered,  the  sub-optimal  control  law  worked  well  for  both  GPS  Block 
HR  and  Hubble  Space  Telescope  in  the  unperturbed  environment.  However,  due  to  the 
inherently  open-loop  nature  of  the  control  law  as  developed,  it  failed  to  perform  well  in 
the  presence  of  perturbing  torques  for  either  vehicle.  The  following  section  presents  some 
specific  ideas  for  future  research  that  will  enhance  the  utility  of  this  control  law  in  real 
applications. 

6.2  Recommendations  for  Further  Study 

The  sub-optimal  control  law  as  currently  developed  suffers  from  two  basic  weaknesses. 
The  first  is  that  vehicle  kinematics  are  not  considered  in  the  formulation,  resulting  in  an 
essentially  arbitrary  (although  small)  rotation  about  the  angular  momentum  vector  at  the 
completion  of  maneuvers,  even  in  the  unperturbed  case.  This  certainly  limits  its  utility 
in  real  applications  where  precise  control  over  pointing  of  communications  and/or  remote 
sensing  hardware  is  required.  Future  work  could  thus  focus  on  including  the  quaternions 
(or  other  attitude  parameters)  into  the  development  of  trajectories  along  the  surface  of 
the  stationary  platform  condition  wheel  momenta  ellipsoid.  The  resulting  torque  control 
trajectories  would  result  in  “smooth”  large  angle  reorientations  with  reduced  attitude 
errors.  Calculation  of  such  trajectories  would  in  all  likelihood  not  be  as  simple  as  the  one 
developed  herein,  and  thus  care  should  be  taken  to  avoid  complicated  computations  that 
could  not  easily  be  performed  by  the  on-board  computer. 

Another  limitation  of  this  control  law  is  that  it  is  based  on  an  assumption  of  perfect 
knowledge  of  the  vehicle  momentum  (x)  and  wheel  momenta  (/i),  and  is  calculated  as¬ 
suming  that  total  angular  momentum  of  the  vehicle  does  not  change.  As  demonstrated  in 
Chapter  5,  external  perturbing  torques  cause  the  angular  momentum  to  change  with  time, 
resulting  in  increased  body  motion  and  final  attitude  and  state  errors.  The  control  law 
thus  needs  to  be  modified  to  allow  for  adjustment  of  control  torques  based  on  actual  vehicle 
angular  momentum.  This  would  require  an  increase  in  the  fidelity  of  the  spacecraft  models 
to  include  vehicle  and  wheel  angular  velocity  estimation  from  sensor  data  (star  sensors. 
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horizon  sensors,  wheel  tachometers,  etc.).  The  attitude  control  logic  could  thus  calculate 
current  angular  momentum,  body  rates,  wheel  angular  momentum,  and  spacecraft  atti¬ 
tude,  constantly  determining  the  wheel  torque  trajectory  necessary  to  move  “closer”  to  the 
stationary  platform  momenta  ellipsoid  (to  regulate  angular  rates)  as  well  as  to  accurately 
reorient  the  vehicle  to  its  final  stationary  attitude. 

Once  a  sub-optimal  control  law  is  found  that  accounts  for  final  attitude,  it  would  be 
interesting  to  compare  the  resulting  trajectory  to  optimal  control  solutions  of  the  problem. 
Specifically,  if  one  solved  an  optimal  control  problem  with  an  objective  function  involving 
the  final  states  (x/,  fij,  and  q/),  as  well  as  weighted  penalties  on  state  variations  and 
control  usage,  it  would  be  interesting  to  see  the  results.  Three  cases  would  be  of  particular 
interest.  First,  specifying  the  end-time  to  be  equal  to  the  sub-optimal  maneuver  would 
allow  direct  comparison  of  state  variations  for  the  two  approaches.  Second,  the  problem 
could  be  formulated  to  obtain  a  minimum  time  solution  to  see  what  the  “fastest”  trajec¬ 
tory  is  (again  compare  state  variations).  Finally,  leaving  end-time  free  would  result  in 
minimization  of  states  and  controls  without  regard  to  time,  and  the  relative  magnitude  of 
the  resulting  angular  velocities  as  compared  to  the  sub-optimal  approach  would  be  very 
interesting.  Hopefully,  results  would  show  that  the  sub-optimal  approach  is  “close”  to  op¬ 
timal  control  solutions  (for  the  fixed  final  time  and  perhaps  minimum  time  cases),  and  by 
virtue  of  its  simplicity,  a  candidate  for  implementation  in  future  attitude  control  systems. 
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Appendix  A.  Constants  and  Unit  Conversions 


Physical  Constants 


3.986012  X  10®  kg  m^fs^ 
4.644  X  10-®  iV/m2 
6378.135  km 


Earth  gravitational  parameter  [3:429] 

Solar  radiation  pressure  at  1  Earth  radius  [1:134] 
Mean  equatorial  radius  of  Earth 


Unit  Conversions 


Source: 

CRC  Tables  [4] 

1  in.  = 

0.0254  m 

1  nmi.  = 

1.15077945  mi. 

1  mi.  = 

1.609344  km 

lib.  = 

0.45359237  kg 

11b.  = 

4.448  N 

A-1 


Appendix  B.  Global  Positioning  System  (GPS)  Block  HR  Modeling  Data 

Sources  for  data  are  indicated  in  the  section  headings.  Data  within  a  section  that  were 
either  derived  or  taken  from  other  sources  are  labeled  accordingly. 

B.l  Mass  and  Inertia  Properties  Calculations 
Source:  [17] 

The  original  mass  properties  data  contained  in  the  GPS  Block  HR  mass  properties 
report  [17]  included  detailed  identification  of  spacecraft  components,  their  locations  within 
the  vehicle,  and  their  masses.  Thus,  foreign  export  of  the  document  was  restricted.  How¬ 
ever,  the  data  extracted  for  use  and  used  in  this  appendix  is  neither  exactly  replicated  nor 
in  sufficient  detail  to  warrant  such  restrictions. 

Figure  B.l  represents  the  fully  deployed  configuration  for  the  GPS  Block  IIR  vehicle 
used  in  mass  properties  calculations.  The  following  data  were  estimated  based  on  source 
data: 


Vehicle  total  mass,  on-orbit  beginning  of  life:  mtot  =  2423.4  lb  ^  1099.2  kg 
Solar  array  mass  (each):  m^a  =  103.304  lb  «  46.858  kg 

Solar  array  surface  area  (each):  Asa  =  7.65  m? 

Solar  array  mass/area  ratio  (derived):  6.125  kgjm? 

Moments  and  products  of  inertia  (in  Tt  frame),  relative  to  spacecraft  center  of  mass, 
on-orbit,  fully  deployed,  beginning  of  life  are  given  as 


hi 

=  15267  in  Ibf  s^  « 

1724.9  kg 

I22 

=  6924  in  Ibf  s^  « 

782.27  kg  m? 

I33 

=  17875  in  Ibf  s"^  w 

2019.5  kg  m? 

I12 

=  hi  = -28S  in  Ibf  s^  « 

-26.324  kg  m? 

Il3 

=  hi  =  —4  in  Ibf  s^  « 

-0.45192  kg  rn 

I23 

=  h2  =  —233  in  Ibf  s^  fa 

—26.32  kg  m? 
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B.2  Material  Properties 

Assumed  properties  for  solar  pressure  modeling  are  as  follows: 


Pa  ^ 
Ps  = 
Pd  = 


0.3  (fraction  of  incident  photons  absorbed) 

0.4  (fraction  of  incident  photons  specularly  reflected) 
0.3  (fraction  of  incident  photons  diffusely  reflected) 


B,3  Momentum  Wheel  Assembly  Properties  Calculations 
Sources:  [19], [7] 

Momentum  wheel  orientation  in  vehicle,  with  respect  to  body  fixed  frame  Ph  are  given  as 


ai 


a2 

as 

a4 


(y2/2  0  - 

(o  -  V2/2 
(-V2/2  0 
(0  V^/2  - 


Figure  B.2  GPS  Block  HR  Wheel  Orientations 
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Axial  wheel  moments  of  inertia 


=  7,2  =  ^3  =  isA  =  0.075  in  lb  ^  0.00847  kg 


Maximum  momentum  wheel  angular  rate:  ±9000  rpm  «  ±943  rad/s 

Minimum  commandable  wheel  motor  torque  (derived):  ±2.768  x  10“'*  N  m 
Maximum  commandable  wheel  motor  torque  (derived):  ±0,071  N  m 


B.4  Orbital  Parameters 
Source:  [18] 

All  of  the  following  parameters  were  assumed  constant  for  modeling  simplicity: 


a 

e 

i 

U) 

n 

T 


26559.91  km  (semi-major  2Lxis) 

0.0039  (eccentricity) 

54.28°  (inclination) 

0.0°  (right  ascension  of  node  (assumed  value  for  modeling  simplicity)) 

7r/2  (argument  of  perigee  (assumed  value  for  modeling  simplicity)) 

0  s  (time  of  perigee  passage  (assumed  value  for  modeling  simplicity)) 


Appendix  C.  Hubble  Space  Telescope  Modeling  Data 

Sources  for  data  are  indicated  in  the  section  headings.  Data  within  a  section  that  were 
either  derived  or  taken  from  other  sources  are  labeled  accordingly. 

C.l  Mass  and  Inertia  Properties  Calculations 
Source:  [2] 

Figure  C.l  represents  the  fully  deployed  configuration  for  the  Hubble  Space  Telescope 
vehicle  used  in  mass  properties  calculations.  The  following  data  were  obtained  from  the 
source  or  assumed  as  indicated: 

Vehicle  total  mass,  on-orbit  beginning  of  life:  mtot  =  24000  lb  «  10886  kg 

Estimated  solar  array  surface  area  (each):  Aga  =  16 

Assumed  solar  array  mass/area  ratio  (from  Appendix  A):  6.125  kgfm^ 

Estimated  solar  array  mass  (derived):  m,a  —  6.1252  (A,a)  ~  98  kg 

The  moments  and  products  of  inertia  were  derived  by  first  calculating  the  inertia 
matrices  for  each  of  the  four  components  of  the  vehicle  (small  cylinder,  large  cylinder,  two 
solar  arrays)  about  their  mass  centers,  calculating  the  mass  center  for  the  entire  vehicle, 
and  then  translating  each  component’s  inertia  matrix  to  the  vehicle  mass  center  using  the 
parallel  axis  theorem  as  outlined  by  Likins  [16:524-525]. 

The  spacecraft  body  was  modeled  as  a  small  right  circular  cylinder  (radius  1.5  m, 
height  4  m)  adjoined  to  a  large  right  circular  cylinder  (radius  2  m,  height  5  m),  both 
of  uniform  mass  distribution.  Each  solar  array  was  modeled  as  a  thin  rectangular  plate 
(height  8  m,  width  2  m,  thickness  0.05  m).  The  solar  arrays  were  assumed  to  be  equal 
distance  from  the  63  axis  and  in  the  b2  —  ba  plane. 
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Given  these  assumptions,  the  moments  and  products  of  inertia  for  the  composite 
spacecraft  about  the  center  of  mass,  expressed  in  Ti,  are 

111  ~  79836  kg  mi 

122  ~  76634  kg  mi 
/33  K,  21679  kg  mi 

112  =  I21  —  0 
7i3  =  I31  =  0 

123  =  I32  =  0 


C.2  Material  Properties 

Assumed  properties  for  solar  pressure  modeling  are  as  follows: 


Pa 

Ps 

Pd 


0.3  (fraction  of  incident  photons  absorbed) 

0.4  (fraction  of  incident  photons  specularly  reflected) 
0.3  (fraction  of  incident  photons  diffusely  reflected) 
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C.3  Momentum  Wheel  Assembly  Properties  Calculations 


Sources:  [2], [20] 

Momentum  wheel  orientation  in  vehicle,  with  respect  to  body  fixed  frame  iFf,  (assumed): 

ai  =  (V2/2  0  -  V2/2Y 
a2  =  (0  -  V2/2  -  V2/2Y 

33  =  {-V2/2  0  -^/2|2Y 

T 

34  =  (0  ^/2/2  -  V^/2) 


Figure  C.2  Hubble  Space  Telescope  Assumed  Wheel  Orientations 


C-4 


Axial  wheel  moments  of  inertia 


isi  =  Is2  =  is3  =  =  0.84  kg  ni^ 


Maximum  momentum  wheel  angular  rate: 

Minimum  commandable  wheel  motor  torque  (Assumed): 
Maximum  commandable  wheel  motor  torque  (Assumed): 


±3000  rpm  «  ±314  rad/s 

0.2  N  m 

3.0  X  10"^  N  m 


C.4  Orbital  Parameters 
Source:  [20] 

All  of  the  following  parameters  were  assumed  constant  for  modeling  simplicity: 


a 

e 

i 

UJ 

n 

T 

o 


6993.135  km 

0.001 

45“ 

0.0” 

7r/2 
0  s 


(semi-major  axis) 

(eccentricity) 

(inclination) 

(right  ascension  of  node  (assumed  value  for  modeling  simplicity)) 
(argument  of  perigee  (assumed  value  for  modeling  simplicity)) 
(time  of  perigee  passage  (assumed  value  for  modeling  simplicity)) 
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Appendix  D.  MATLAB  Code  Listings 

The  following  programs  were  developed  by  the  author  from  January-November  1995. 
They  are  MATLAB^^  script  files,  and  along  with  the  standard  MATLAB^^  library  rou¬ 
tines,  provide  full  functionality  for  all  simulations  conducted  in  this  report.  The  programs 
were  run  on  a  Sun  SPARCstation  20  workstation  under  Sun  Operating  System  ver.  4.1  or 
later  and  MATLAB^^  ver.  4.2c. 

Figure  D.l  portrays  the  functional  flow  of  data  within  the  programs. 


Figure  D.l  Simulation  Program  Functional  Flow  Diagram 
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D.l  GPS  Block  HR  Data  Definition  File 


1 

2 

3 

4 

5 

6 

7 

8 
9 

10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 
21 
22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

36 

37 

38 

39 

40 

41 

42 

43 

44 


y. 

satmod=l;  */,  Satellite  model  to  use  for  3-D  plots  (1=GPS,  2=Hubble) 

pltflgl=l;  */•  Turn  momentum  ellipsoid  plots  on/off  (1/0) 

ellflg=l;  7.  Select  momentum  ellipsoid  plotting  method  and  resolution 
inc0=.01;  7,  0=achievable  regions  of  ellipsoid,  point-by-point  (slow) 

incl=3;  7,  l=entire  ellipsoid,  peorameterized  method  (fast) 

pltflg2=l;  7,  Turn  dimen.  momenta,  ang.  veloc.  plots  on/off  (1/0) 

pltflg3=0;  7,  Turn  non-dimen.  momenta,  ang.  veloc.  plots  on/off  (1/0) 

pltflg4=l;  7,  Turn  3-D  attitude  plots  on/off  (1/0) 

pertflg=0;  7.  Turn  grav.  grad,  and  solar  press,  torques  on/off  (1/0) 

7.  0  =none  1  =sol.  press,  only  2  =sol.  press,  and  grav. 

small=0.1;  7.  Value  of  small  parameter  on  torques.  Should  be  (+) 

r=l;  7,  Momentum  wheels  to  use 

7.  0  =  (1/2/3),  1  =  (1/2/4),  2  =  (1/3/4),  3  =  (2/3/4) 

ph_o=0;  7.  Initial  Euler  angle  phi  about  e3  axis  (deg) 

th_o=0;  7,  Initial  Euler  angle  theta  about  el  axis  (deg) 

ps_o=0;  7,  Initial  Euler  angle  psi  about  e3  eucis  (deg) 

7. 

ph_f=45;  7i  Final  Euler  angle  phi  about  e3  axis  (deg) 

th_f=45;  7.  Final  Euler  angle  theta  about  el  axis  (deg) 

ps_f=0;  7.  Final  Euler  angle  psi  about  e3  axis  (deg) 

7. 

h_o=5;  7.  Magnitude  of  vehicle  ang.  momentum  (constant)  [kg*m"2/s] 

7,  used  to  dimensionalize  results.  Note  that  since  vehicle 
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45 

46 

47 

48 

49 

50 

51 

52 

53 

54 

55 

56 

57 

58 

59 

60 
61 
62 

63 

64 

65 

66 

67 

68 

69 

70 

71 

72 

73 

74 

75 

76 

77 

78 

79 

80 
81 
82 

83 

84 

85 

86 

87 

88 
89 


*/.  ang.  veloc.  "0,  larger  h_0  < — >  higher  rotor  momenta; 
mu_o=[l  0  0];  */,  Initial  non-dimen.  rotor  momenta  (stationary  platform) 

*/,  format :  mu_o=  [mul  mu2  mu3] 

7,  use  mu_o=-99  to  generate  random  point 
ths0=135;  7,  Initial  sun  vector  angle  from  el  (CCW,  in  degrees) 

traj=2;  7.  Wheel  torque  trajectory 

7.  (l=direct  2=suboptimal) 
t0=0;  7.  Start  time  for  integration  (sec). 

tf=-999;  7t  Stop  time  for  integration  (sec).  Use  -999  for  autocalc. 

7.  NOTE:  non-dimen.  values  calculated  below  s/c  properties 
num=2;  7  desired  number  of  3-D  attitude  plots  for  this  time  span 

7. 

7. 


7.77.7.7.  SPACECRAFT  PHYSICAL  PROPERTIES  7.77.7.7. 


7.  Inertia,  Mass,  and  coordinate  frame  data  obtained  from  Navstar  GPS 
7.  Block  HR  Phase  II  Mass  Properties  Report,  9  Mar  95 

7. 


al=[sqrt(2)/2;  0; 

a2=[0;  -sqrt(2)/2; 

a3=[-sqrt(2)/2;  0; 

a4=[0;  sqrt(2)/2; 

Isl=. 00847; 

Is2=. 00847; 

Is3=. 00847; 

Is4=.00847; 
ws lmax=9000*2*pi/ 60 ; 
ws2max=9000*2*pi/60 ; 
ws3max=9000*2*pi/60 ; 
ws4max=9000*2*pi/60 ; 

7. 

I=[  15267  -233  -4; 

-233  6924  -46; 

-4  -46  17875] ; 


-sqrt(2)/2];  7.  wheel  orientation  w.r.t 

-sqrt(2)/2];  7  body  fixed  frame  b  (and  XYZ) 

-sqrt(2)/2] ; 

-sqrt(2)/2] ; 

7  Axial  moments  of  inertia  of  the 
7  4  momentum  wheels  (kg*m~2) 


7  Maximum  wheel  angular  speed  magnitudes  (rad/s) 


GPS  Blk  HR  on-orbit  inertia  tensor, 
relative  to  mass  center  in  XYZ  frame 
units  are  in*lb*sec“2 


1= (.0254)* (4. 448) *1; 


Convert  from  in*lb*sec“2  to  kg*m“2 
(.0254  m/ in) (4.448  N/lb) (in*lb*sec"2)= 
N*m*sec“2  =  (kg*m/sec~2)*m*sec''2=kg*m'‘2 


Ic=trace(I) ; 

I_d=I ; 
I=I_d./Ic; 


7.  Charact.  inertia  for  non-dimensionalizing 
7.  (note:  this  is  frame  invariant) 

7  Non  dimens ionalized  inertia  matrix 
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90 

91 

92 

93 

94 

95 

96 

97 

98 

99 
100 
101 
102 

103 

104 

105 

106 

107 

108 

109 

110 
111 
112 

113 

114 

115 

116 

117 

118 

119 

120 
121 
122 

123 

124 

125 

126 

127 

128 

129 

130 

131 

132 

133 

134 


•/. 

Is_d=eye(3) ; 
Is_d(l,l)=Isl; 
Is_d(2,2)=Is2; 
Is_d(3,3)=Is3; 
if  r==l 

Is_d(3,3)=Is4; 
elseif  r==2 
Is_d(2,2)=Is3; 
Is_d(3,3)=Is4; 
elseif  r==3 
Is_d(l,l)=Is2; 
Is_d(2,2)=Is3; 
Is_d(3,3)=Is4; 
end; 

•/. 

Is=Is_d./Ic; 

•/. 


'/,  Create  diagonal  matrix  template; 
*/,  Default  Configuration  (r=0) ; 


y,  Fill  remaining  diagonal  entries  of  Is 
*/,  Based  on  choice  of  wheel  combination. 


y,  Non  dimens ionalized  inertia  of  wheels 


y,  Now  define  vectors  from  vehicle  c.m.  to  component  centers 

y,  of  pressure  for  solar  pressure  modeling  (m)  ,  the  effective  areas  (m~2) 

y,  and  the  surface  normal  vectors  (all  in  Fb) 


y. 

rs(l, :)=[!  0  0] ; 
rs(2, : )  =  [-!  0  0]  ; 
rs(3,:)=[0  1  0]; 
rs(4,:)  =  [0  -1  0]; 
rs(5, :  )  =  [0  0  -1]  ; 
rs(6, :  )  =  [0  0  1]  ; 
rs(7, :)=[0  3.5  .37] ; 
rs (8 , : ) = [0  -3 . 5  . 37] ; 
rs(9, :)=[0  3.5  .37] ; 
rs(10, :)=[0,-3.5  .37]; 
As(l)=4.424; 

As (2) =4. 424; 

As (3) =4. 424; 

As (4) =4. 424; 

As (5) =2. 4964; 

As (6) =2. 4964; 
As(7)=7.65; 

As (8) =7. 65; 

As(9)=7.65; 
As(10)=7.65; 
ns(l, : )  =  [1  0  0]  ; 
ns(2, :)  =  [-!  0  0]  ; 
ns(3,:)=[0  1  0]; 


y.  +bl 

y.  -bi 
y.  +b2 
y.  -b2 
y.  -b3 
y.  +b3 
y.  +b2 
y.  -b2 
y.  +b2 
y.  -b2 


face 

face 

face 

face 

face 

face 

solar 

solar 

solar 

solar 


array  FRONT 
array  FRONT 
array  BACK 
array  BACK 
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135 

ns  (4, 

)=[0 

- 

0]; 

136 

ns  (5, 

O 

i_j 

11 

0 

-1]; 

137 

ns  (6, 

II 

1 — 1 
o 

0 

1]; 

138 

ns  (7, 

)=[1 

0 

0]; 

139 

ns  (8, 

V— 1 

1 _ 1 

II 

0 

0]; 

140 

ns  (9, 

)=[-. 

L  ( 

0  0]; 

141 

ns (10, 

,:)  =  [■ 

“1 

0  0]; 

142 

y. 

143 

144  VLW.  ORBITAL  PARAMETERS 

145 

146  '/. 

147  '/,  Note:  right  ascension  of  the  node  is  assumed  to  be  0, 

148  '/,  and  the  argument  of  perigee  is  assumed  to  be  pi/2 

149  y. 

150  sma=26559.91;  '/,  Semi-major  axis  (km) 

151  ecc=0.0039;  '/,  Eccentricity 

152  incl=deg2rad(54.28) ;  '/.  Inclination  (degrees) 

153  y. 

154 

155  y.y.y.y.y.  non-dimensionalize  times  y.y.y.y.y. 

156 

157  y. 

158  y,  This  is  necessary  for  equation  consistency  for  integration 

159  y.  Must  be  done  here  since  Ic  needed  from  above 

160  y. 


y,  Note:  right  ascension  of  the  node  is  assumed  to  be  0, 
y,  and  the  argument  of  perigee  is  assumed  to  be  pi/2 

y. 

sma=26559.91;  '/,  Semi-major  axis  (km) 

ecc=0.0039;  '/,  Eccentricity 

incl=deg2rad(54.28)  ;  '/.  Inclination  (degrees) 

y. 

•/•/•/•/•/•/•/•/ V  •/•/•/ V  VV  •/•/•/•/ VV  VV  VV  VV  •/ V  •/•/•/•/•/•/•/ V  •/ V  V  •/•/•/ V  •/ V  vv 


y.y.y.y.y.  non-dimensionalize  times  y.y.y.y.y. 


161  t0=h_o*t0/lc; 

162  if  tf>0,  tf=h_o*tf/Ic;  end;  */,  don’t  recalculate  if  tf=-999 
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D.2  Hubble  Space  Telescope  Data  Definition  File 


1 

2 

3 

4 

5 

6 

7 

8 
9 

10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 
21 
22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

36 

37 

38 

39 

40 

41 

42 

43 

44 


•/. 

satmod=2;  */,  Satellite  model  to  use  for  3-D  plots  (1=GPS,  2=Hubble) 

pltflgl=0;  */,  Turn  momentvim  ellipsoid  plots  on/off  (1/0) 

ellflg=l;  */,  Select  momentum  ellipsoid  plotting  method  and  resolution 
inc0=.05;  ’/,  0=achievable  regions  of  ellipsoid,  point-by-point  (slow) 

incl=3;  '/,  l=entire  ellipsoid,  pcirameterized  method  (fast) 

pltflg2=0;  y,  Turn  dimen.  momenta,  ang.  veloc.  plots  on/off  (1/0) 

pltflg3=0;  y,  Turn  non-dimen.  momenta,  ang.  vel.  plots  on/off  (1/0) 

pltflg4=l;  y,  Turn  3-D  attitude  plots  on/off  (1/0) 

pertflg=l;  '/,  Turn  grav.  grad,  and  sol.  press,  torques  on/off 

y,  0  =none  1  =sol.  press,  only  2  =sol.  press,  and  grav. 
small=.5;  */,  Value  of  small  parameter  on  torques.  Should  be  (+) 

r=0;  y.  Momentum  wheel  triads  to  use 

y.  0  =  (1/2/3),  1  =  (1/2/4),  2  =  (1/3/4),  3  =  (2/3/4) 
ph_o=0;  y.  Initial  Euler  angle  phi  about  e3  axis  (deg) 

th_o=0;  y,  Initial  Euler  angle  theta  about  el  axis  (deg) 

ps_o=0;  y.  Initial  Euler  angle  psi  about  e3  axis  (deg) 

y. 

ph_f=0;  y.  Final  Euler  eingle  phi  about  e3  axis  (deg) 

th_f=90;  y.  Final  Euler  angle  theta  about  el  axis  (deg) 

ps_f=0;  y,  Final  Euler  angle  psi  about  e3  axis  (deg) 

y. 

h_o=200;  y.  Magnitude  of  vehicle  ang.  momentum  (constant)  [kg*m~2/s] 

y,  used  to  dimensionalize  results.  Note  that  since  vehicle 


45 

46 

47 

48 

49 

50 

51 

52 

53 

54 

55 

56 

57 

58 

59 

60 
61 
62 

63 

64 

65 

66 

67 

68 

69 

70 

71 

72 

73 

74 

75 

76 

77 

78 

79 

80 
81 
82 

83 

84 

85 

86 

87 

88 
89 


*/,  ang.  vel.  are  '0,  larger  h_0  < — >  higher  rotor  momenta; 
mu_o=Cl  0  0];  */,  Initial  non-dimen.  rotor  momenta  (stationary  platform) 

'/,  format:  mu_o=[mul  mu2  mu3] 

'/,  use  mu_o=-99  to  generate  random  point 
ths0=135;  */,  Initial  sun  vector  angle  from  el  (CCW,  in  degrees) 

traj=2;  */,  Wheel  torque  trajectory 

*/,  (l=direct  2=suboptimal) 

t0=0;  */.  Start  time  for  integration  (sec).  Use  -999  for  autocalc. 

tf=-999;  y.  Stop  time  for  integration  (sec).  Use  -999  for  autocalc. 

y,  NOTE:  non-dimen.  values  calc,  below  s/c  properties 
num=20;  '/,  desired  number  of  3-D  attitude  plots  for  this  time  span 

y. 

y. 


y.y,y.y.y.  spacecraft  physical  properties  y.y.y.y.y, 


y. 

y,  The  Hubble  Space  Telescope  is  modeled  as  two  right  circular  cylinders 
y,  of  uniform  density,  with  two  rectangular  parallelepiped 
y,  solar  arrays  at  equal  disteinces  from  s.c.  c.m. 

y. 

y,  there  are  four  axisymmetric  momentum  wheels 

y. 

al=[sqrt(2)/2;  0;  -sqrt(2)/23 ;  '/,  wheel  orientation  w.r.t 

a2=[0;  -sqrt(2)/2;  -sqrt(2)/2];  '/,  body  fixed  frame  b  ) 

a3=[-sqrt(2)/2;  0;  -sqrt(2)/2]; 

a4=[0;  sqrt(2)/2;  -sqrt(2)/2]; 


Isl=.84 


y,  Axial  moments  of  inertia  of  the 


Is2=.84 


y,  4  momentum  wheels  (kg*m~2) 


Is3=.84 


Is4=.84; 

ws lmax=3000*2*pi/60 ; 

ws2max=3000*2*pi/60 ; 

ws3max=3000*2*pi/60 ; 

ws4max=3000*2*pi/ 60 ; 

m_tot=24000/2 . 204622622 ; 

r_sm=1.5; 

h_sm=4; 

r.lg=2; 

h_ig=5: 

h_sa=8; 

w_sa=2; 

d_sa=.05; 

rho_d=.3; 

rho_s=.4; 


y,  Max.  wheel  ang.  speed  magnitudes  (rad/s) 


y,  HST  total  mass  (lb) 

y,  Radius  of  small  cylinder  (m) 

y.  Height  of  small  cylinder  (m) 

y,  Radius  of  large  cylinder  (m) 

y.  Height  of  large  cylinder  (m) 

y,  Height  of  each  solar  array  (m) 

y,  Width  of  each  solar  eirray  (m) 

y,  Thickness  of  each  solar  array  (m) 

y.  Spacecraft  surface  diffuse  reflectivity 

y  Spacecraft  surface  specular  reflectivity 


90  •/. 

91  gps_sa_in2a=6 . 125;  */. 

92  m_sa=h_sa*w_sa*gps_sa_m2a;  */, 

93  m_sc=m_tot-2*m_sa;  7, 

94  vol_sm=pi*r_sm*2*h_sm;  7, 

95  vol_lg=pi*r_lg“2*h_lg;  7. 

96  f  rac_sm=vol_sm/  (vol_sm+vol_lg)  ;  7. 

97  m_sm=frac_sm*m_sc;  7. 

98  m_lg=m_sc-m_sm;  7# 

99  7. 

100  d=[in_siii  m_lg  2*m_sa:  1  -1  0;  0  -1 

101  R_sm=[0  0  d(l)];  7. 

102  R_lg=[0  0  d(2)];  7. 

103  R_sl=[0  4  d(3)]  ;  7. 

104  R_s2=[0  -4  d(3)] ; 

105  7. 

106  I_sm=eye(3);  7. 

107  I_sm(l , l)  =  (3*r_sm~2+h_sm*2)/12;  7. 

108  I_sm(2,2)=I_sm(l,l);  */. 

109  I_sm(3,3)=r_sm“2/2; 

110  I_sm=in_sm*I_sm; 

111  7. 

112  I_lg=eye(3);  7. 

113  I_lg(l,l)  =  (3*r_lg-2+h_lg~2)/12;  7. 

114  I_lg(2.2)=I_lg(l,l);  7. 

115  I_lg(3,3)=r_lg“2/2; 

116  I_lg=m_lg*I_lg; 

117  7. 

118  I_sa=Gye(3); 

119  I_sa(l,l)  =  (h_sa"2+w_sa'‘2)/12;  7. 

120  I_sa(2,2)  =  (h_sa~2+d_sa“2)/12;  7. 

121  I_sa(3,3)  =  (d_sa“2+w_sa“2)/12;  7. 

122  I_sa=in_sa*I_sa: 

123  7. 

124  R_smx=[  0  -R_sm(3)  R_sm(2) ;  7. 

125  R_sm(3)  0  -R_sm(l);  % 

126  -R_sm(2)  R_sm(l)  0]; 

127  R_lgx=[  0  -R_lg(3)  R_lg(2);  7. 

128  R_lg(3)  0  -R.lgd);  7. 

129  -R_lg(2)  R_lg(l)  0]; 

130  R_slx=[  0  -R_sl(3)  R_sl(2);  7. 

131  R_sl(3)  0  -R.sKl);  7. 

132  -R_sl(2)  R.sKl)  0]; 

133  R_s2x=[  0  -R_s2(3)  R_s2(2);  7. 

134  R_s2(3)  0  -R_s2(l);  7. 


GPS  HR  solar  array  mass/area  ratio 
Mass  of  each  HST  solar  array  (kg) 

Mass  of  HST  without  solar  arrays  (kg) 
Volume  of  small  cylinder  (m"3) 

Volume  of  large  cylinder  (m“3) 

Fraction  of  total  voliuae  from  small  cyl. 
Mass  of  small  cylinder  (kg) 

Mass  of  large  cylinder  (kg) 

1]\[0; (h_sm+h_lg)/2;h_lg/2] ; 

Vectors  (in  b  frame)  from  system  c.m, 
to  small  and  large  cylinder  c.m.  (m) 
and  solar  array  c.m. 

Inertia  matrix  for  small  cylinder  about 
center  of  mass,  expressed  in  b  frame 
(kg*m*2) 


Inertia  matrix  for  large  cylinder  about 
center  of  mass,  expressed  in  b  frame 
(kg*m''2) 


Inertia  matrix  for  solar  array  about 
center  of  mass,  expressed  in  b  frame 
(kg*m'‘2) 

Skew  S3rinmetrix  matrix  for  translating 
I_sm  to  system  c.m. 

Skew  symmetrix  matrix  for  translating 
I.lg  to  system  c.m. 

Skew  symmetrix  matrix  for  translating 
I^sal  to  system  c.m. 

Skew  symmetrix  matrix  for  translating 
I.sa2  to  system  c.m. 
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135 

-R_s2(2)  R.s2(l)  0]; 

136 

I _  sm_  c=I _  sm-m_  sm*R_  smx*R_  smx ; 

y.  Translate  inertia  matricesto  system 

137 

I_lg_c=I_lg-m_lg*R_lgx*R_lgx ; 

7  center  of  mass,  in  b  frame 

138 

I_sal_c=I_sa-m_sa*R_slx*R_slx; 

139 

I_sa2_c=I.sa-m_sa*R^s2x*R_s2x; 

140 

I=I_sm.c+I_lg_c+I_sal_c+I.sa2_c;  7,  Construct  total  inertia  matrix  about 

141 

y,  c.m.  in  b  frame 

142 

Ic=trace(I) ; 

7  Charac.  inertia  for  non-dimensionalizing 

143 

7  (note:  this  is  frame  invariant) 

144 

I_d=I ; 

y,  dimensional  inertia  matrix 

145 

1=1. /Ic; 

y.  Non  dimens ionalized  inertia  matrix 

146 

y. 

147 

Is_d=eye(3) ; 

y,  Create  diagonal  matrix  template; 

148 

Is_d(l,l)=Isl; 

y,  Default  Configuration  (r=0) ; 

149 

Is_d(2,2)=Is2; 

150 

Is_d(3,3)=Is3; 

151 

if  r==l 

y,  Fill  remaining  diagonal  entries  of  Is 

152 

Is_d(3,3)=Is4; 

y,  Based  on  choice  of  wheel  combination. 

153 

elseif  r==2 

154 

Is_d(2,2)=Is3; 

155 

Is_d(3,3)=Is4; 

156 

elseif  r==3 

157 

Is_d(l,l)=Is2; 

158 

Is_d(2,2)=Is3; 

159 

Is_d(3,3)=Is4; 

160 

end; 

161 

y. 

162 

Is=Is_d./Ic; 

7  Non  dimensionalized  inertia  of  wheels 

163 

y. 

164 

7  Now  define  vectors  from  vehicle  c.m.  to  component  centers 

165 

7  of  pressure  for  solar  pressure  modeling  (m) ,  the  effective  areas  (m*2) 

166 

y,  and  the  surface  normal  vectors  (all  in  Fb) 

167 

y. 

168 

rs(l,:)  =  [0  0  0];  '/,  small 

cylinder.  Dummy  value,  calc,  during  integr. 

169 

rs(2,:)  =  [0  0  0];  */,  large 

cylinder.  Dummy  value,  calc,  during  integr. 

170 

rs(3, :)=R_sm+[0  0  h_sm/2] ; 

7  small  cylinder  top 

171 

rs(4,:)=R_lg+[0  0  -h_lg/2] ; 

y,  laurge  cylinder  bottom 

172 

rs(5, ;)=R_sl; 

y,  +b2  solar  array  FRONT 

173 

rs(6, :)=R_s2; 

y,  -b2  solar  array  FRONT 

174 

rs(7, :)=R_sl; 

y,  +b2  solar  array  BACK 

175 

rs(8, :)=R.s2; 

y,  -b2  solar  array  BACK 

176 

rs(9, :)  =  [0  0  0]  ; 

7  dummy  value,  no  9th  surface 

177 

rs(10, :)=[0  0  0] ; 

y,  dummy  value,  no  10th  surface 

178 

As ( 1 ) =2*r_sm*h_sm ; 

179 

As (2) =2*r_lg*h_lg ; 
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180  As(3)=pi*r_sm“2; 

181  As(4)=pi*r_lg“2; 

182  As(5)=h_sa*w_sa; 

183  As(6)=h_sa*w_sa; 

184  As(7)=h_sa*w_sa; 

185  As(8)=h_sa*w_sa; 

186  As(9)=0; 

187  As(10)=0; 


188 

ns(l. 

)  =  [0 

0  0]; 

189 

ns  (2, 

)  =  [0 

0  0]; 

190 

ns  (3, 

)  =  [0 

0  1]; 

191 

ns  (4, 

)=[o 

0  -1]; 

192 

ns  (5, 

)  =  [1 

0  0]; 

193 

ns  (6, 

)  =  C1 

0  0]  ; 

194 

ns  (7, 

)  =  [-: 

10  0]; 

195 

ns  (8, 

)  =  [-: 

10  0]; 

196 

ns  (9, 

)  =  [0 

0  0]; 

197 

ns (10, 

,:)  =  [0  0  0]; 

'/,  small  cylinder.  Dummy  value,  calc,  during  integration 
y,  large  cylinder.  Dummy  value,  calc,  during  integration 


198 

199 

200 
201 
202 

203 

204 

205 

206 

207 

208 

209 

210 
211 
212 

213 

214 

215 

216 

217 

218 

219 

220 


y. 

y. 


y.y.y.y.y.  orbital  parameters  y.y.y.y.y. 

y. 

y.  Note:  right  ascension  of  the  node  is  assumed  to  be  0, 
y,  and  the  argument  of  perigee  is  assumed  to  be  pi/2 

y. 

sma=6993 . 135 :  */  Semi-major  axis  (km) 

ecc=0.001;  */.  Eccentricity 

incl=deg2rad(45)  ;  '/,  Inclination  (degrees) 

y. 

y. 


NON-DIMENSIONALIZE  TIMES 


Vim 


y.y.y.y.y. 

vimmmmmi:immm:imm:imm 

y. 

y,  This  is  necessary  for  equation  consistency  for  integration 
y,  Must  be  done  here  since  Ic  needed  from  above 

y. 

t0=h_o*t0/lc; 

if  tf>0,  tf=h_o*tf /Ic;  end;  */,  don’t  recalculate  if  tf=-999 
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D.3  Main  Simulation  Program 


rotor.m.m 

Capt  Greg  Schultz,  GA-95D 


•/.•///. 

va 
m 


revision 


11.02.95.1230 


1 

2  •///.•/. 

3  •/.*/.•/. 

4  •///.•/. 

5  •/.•///. 

6  •/.*/.•/. 

7  •/,•/.•/, 

8 

9  •/,•/.*/. 

10  m 

11 

12  y. 

13  clear; 

14  format  short  e; 

15  global  I_d  Ic  Is_d  h_o  JI  small; 

16  global  A  B  alal  ala2  a2a2; 

17  global  mu_o  mu_f  Rgmu  Rdmu; 

18  global  r  direct  wsmax  ellflg  GM; 

19  global  rho_d  rho_s; 

20  global  rs  ns  As  r_sm  r_lg  R_sm  R_lg  satmod  pertflg; 

21  •/. 

22  File=input(’Data  file  (in  single  quotes  w/o  .m  extension,  e.g.  dat.def):’); 

23  •/. 

24  eval(File);  '/.  Execute  file  to  load  data  for  simulation 

25  y. 

26  y,attpos=  [0  25  1200  900]  ; 

27  attpos=[0  100  700  700]; 

28  y. 

29  if  r==0 
A=  [al  a2  a3]  ; 

wsmax=  [ws  Imaix  ws2m2ix  ws3max]  ; 


y,  display  format 

y,  Define  global  variables  used 

y,  here  and  in  subroutines 


30 

31 

32  elseif  r==l 

33  A=  [al  a2  a4]  ; 

34  wsmax=[wslm3ix  ws2max  ws4max]  ; 

35  elseif  r==2 

36  A=  [al  a3  a4]  ; 

37  wsmax=  [ws  Imeuc  ws3max  ws4mauc]  ; 

38  elseif  r==3 

39  A=[a2  a3  a4] ; 

40  wsmax=  [ws2mEix  ws3max  ws4max]  ; 

41  end; 

42  al=A(:,l);  a2=A(:,2);  a3=A(:,3); 

43  y. 

44  if  rank (A) <3 


y,  Use  this  for  presentations 
y,  [Xo  Yo  Xwidth  Ywidth]  pixels 

y.  Calculate  A  and  wsmax  based  on  wheels 
y,  selected  to  use  in  simulation 


y,  redefine  al..a3  for  subsequent  calcs, 
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45 

46 

47 

48 

49 

50 

51 

52 

53 

54 

55 

56 

57 

58 

59 

60 
61 
62 

63 

64 

65 

66 

67 

68 

69 

70 

71 

72 

73 

74 

75 

76 

77 

78 

79 

80 
81 
82 

83 

84 

85 

86 

87 

88 
89 


’  Wheel  unit  vectors  do  not  form  a  basis  -  modify  data  input  file’ 
break; 
end; 

•/. 


SPECIFY  IMITIAL  AND  FINAL  CONDITIONS 


•/. 

*/,  Convert  initial,  final  Euler  angles  to  radians 

7. 

ph_o=deg2rad(ph_o) ;  th_o=deg2rad(th_o) ;  ps_o=deg2rad(ps_o) ; 
ph_f=deg2rad(ph_f ) ;  th_f =deg2rad(th_f ) ;  ps_f=deg2rad(ps_f ) ; 

7. 

%  Obtain  rotation  matrices  Ro,  Rf  relating  initial,  desired  final 
7,  position  of  body  frame  to  inertial. 

7. 

Ro=rot3(ps_o)*rotl(th_o)*rot3(ph_o) ;  7.  xo|b=Ro*xo|i 
Rf=rot3(ps_f  )*rotl(th_f  )*rot3(ph_f )  ;  7#  xf|b=Rf*xf|i 

7. 

7i  calculate  initial  Euler  axis  (lambda)  and  rotation  angle  (m) 

7. 

[eigvec,eigval]=eig(Ro) ; 
if  abs(l-eigval(l, !))<=. 001 
lambda=eigvec( : , 1) ; 
m=angle(eigval(2,2)) ; 
elseif  abs(l-eigval(2,2))<=.001 
lambda=eigvec( : ,2) ; 
m=angle(eigval(l,l)) ; 
elseif  abs(l-eigval(3,3))<=.001 
lambda=eigvec( : ,3) ; 
m=angle(eigval(l,l)) ; 
else 


errorC ’CANNOT  FIND  EULER  AXIS 
end; 

m=abs (m) ; 

7. 

ql0=lambda(l)*sin(m/2) ; 
q20=lambda(2)*sin(m/2) ; 
q30=lambda(3)*sin(m/2) ; 
q40=cos(m/2) ; 
q0=  CqlO  q20  q30  q40] ; 

7. 

if  mu_o==-99 

xx=rand(size(l :3)) ; 
x_o=xx/norm(xx) ; 


No  eigenvalue  close  to  1’) 

7.  take  pos.  rotations  about  euler  axis 
7.  calculate  initial  quaternions 


7.  generate  random  mu_o  if  desired 
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90 

91 

92 

93 

94 

95 

96 

97 

98 

99 
100 
101 
102 

103 

104 

105 

106 

107 

108 

109 

110 
111 
112 

113 

114 

115 

116 

117 

118 

119 

120 
121 
122 

123 

124 

125 

126 

127 

128 

129 

130 

131 

132 

133 

134 


mu_o=inv(A)*x_o’ ; 
end; 

•/. 

if  sizG(mu_o,l)==l,  mu_o=mu_o’;  end; 

•/. 

mu_f=inv(A)*Rf*Ro’*A*mu_o;  '/,  calc.  mu_f  based  on  init.,  desired 

'/,  final  orientations 

•/. 

x_o=A*mu_o ; 
x_f =A*mu_f ; 

’Dimensionless  s/c  ang.  mom.  (stationary  platform)  boundary  conditions  are:’ 


x_o 

x_f 

’Dimensionless  wheel  ang.  mom.  (stationary  platform)  boundary  conditions  are 

mu_o 

mu_f 

’Initial  wheel  angular  velocities  (rad/sec):’ 
wso=inv(Is*Ic)*h_o*mu_o 

•/. 

J=I-A*Is*A’;  JI=inv(J) ; 

y. 

y,  Calculate  initial  orbital  position,  velocity 
y,  (assuming  right  ascension=0,  arg  per=pi/2) 

y. 

if  pertflg>0 

GM=398601.2;  '/.  Gravitational  parameter 

rper=sma*(l-ecc) ;  */.  Radius  of  perigee 

enr=-GM/(2*sma) ;  '/.  Orbital  energy 

Rorb0=[0  rper*cos(incl)  rper*sin(incl)]  ;  '/.  Initial  orbital  position 

VorbO=[-(2*(enr+GM/rper))“.5  0  0];  '/,  Initial  orbital  velocity 

TP=2*pi*(sma''l  .5)/(GM''  .5)  ; 
end; 

y. 

ths0=deg2rad(ths0)  ;  */,  convert  initial  sun  angle  to  rad 

y. 

y.y.y.y.y.  calculate  direct  wheel  torque  control  y,y,y.y.y. 

y. 

if  traj==l 

direct=small* ( (mu_f -mu_o) / (norm(mu_f -mu_o) ) ) ; 


if  tf==-999 

tf=(norm(mu_f-mu_o)) /small;  '/,  auto-calculate  final  time  for 

end;  '/.  desired  final  condition  mu_f 

end; 
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135  y. 

136 


137  '/;/.•/.  CALCULATE  SUBOPTIMAL  WHEEL  TORQUE  CONTROL 

138 

139  y. 

140  y,  define  new  gi  frame  in  which  origin,  mu_o,  mu_f  lie  in  gl-g2  plane 

141  y,  this  is  a  suboptimal  trajectory  NEAR  to  the  "great  circle" 

142  y,  trajectory  on  a  sphere.  Approximation  is  best  when  ellipsoid  is 

143  y,  nearly  spherical  (rotors  nearly  orthogonal)  . 

144  y,  Only  limitation  is  must  have  non-zero  cross  prod,  between  mu_o  auid  mu_f 

145  y. 

146  if  traj==2 

147  gl=mu_o/norm(mu_o) ; 

148  g3=cross(mu_o,mu_f)/(norm(cross(mu_o,mu_f))) ; 

149  g2=cross(g3,gl)/norm(cross(g3,gl)) ; 

150  Rgmu=[gl  g2  g3] ;  */,  rot.  matrix  from  g  to  "mu"  frame 

151  AG=A*Rgmu; 

152  agl=AG(:,l);  ag2=AG(:,2);  ag3=AG(:,3); 

153  alal=agl’*agl;  a2a2=ag2’*ag2;  ala2=agl’*ag2; 

154  B=[ala2  a2a2  0;-alal  -ala2  0;0  0  0];  '/,  Wheel  torque  in  "nu"  frame 

155  y. 

156  y,  auto-calculate  final  time  for  desired  final  condition  mu_f 

157  y. 

158  if  tf==-999 

Bl=[ala2  a2a2;-alal  -ala2]  ;  '/.  Two  rotor  torque  law  in  gi  frame 

nu_f =Rgmu ’  *mu_f ;  '/.  rotor  momenta  in  gi  frame 

nu_f=nu_f  (1:2) ;  '/.  in  gi  frame,  nu_f(3)=0  anyway 

[vl,dl]=eig(Bl) ; 

Pbl=[real(vl(: ,1))  imag(vl(: ,1))] ; 
wl=imag(dl(l,l)) ; 
abl=Pbl\Cl;0] ; 
aal=abl(l);  bl=abl(2); 
vvl=Pbl\nu_f ; 
vll=vvl(l);  vl2=vvl(2); 

tf=atan2(bl*vll-aal*vl2,  aal*vll+bl*vl2)/(wl*small) ; 
if  tf<0,  tf=-tf;  end; 


159 

160 
161 
162 

163 

164 

165 

166 

167 

168 

169 

170 

171  end; 

172  end; 

173  y. 


174 

175  y.y.y.y.y.  integrate  equations  of  motion  y,y,y.y,y, 

176  nnvixixixm^^ 

177  y. 

178  t0_d=lc*t0/h_o;  '/,  dimens ionalize  integration  times 

179  tf_d=Ic*tf/h_o: 


D-14 


180  ’integration  time  (sec):’ 

181  t0_d 

182  tf_d 

183  •/. 

184  if  pertflg>01 

185  ’Number  of  orbital  periods  spanned  by  this  maneuver  :  ’ 

186  tf_d/TP 

187  end; 

188  */,  unperturbed  RHS  file,  tolerances,  i.c.s 

189  if  pertflg==0 

190  if  traj==l,  File=’rotor_el’ ; 

191  elseif  traj==2,  File=’rotor_e2’ ; 

192  end ; 

193  tol=le-10 

194  statesO=Cx_o’  mu_o’  qO]  ; 

195  [t, states] =ode45 (File, to, tf,states0,tol,0) ; 

196  end ; 

197  */,  perturbed  RHS  file,  tolerances,  i.c.s 

198  if  pertflg>0 

199  if  traj==l,  File=’rotor_el_per’ ; 

200  elseif  traj==2,  File=’rotor_e2_per’ ; 

201  end; 

202  tol=le-10 

203  states0=[x_o’  mu_o’  qO  thsO  RorbO  VorbO] ; 

204  [t, states] =ode45 (File, to, tf, St ates0,tol,0) ; 

205  end; 

206  •/. 

207  */,  integration  is  now  done,  so  let’s  redefine  returned 

208  */,  state  variables  with  more  recognizable  names 

209  •/. 

210  */,  These  are  returned  for  both  perturbed  and  unperturbed  system 

211  •/. 

212  xl  =states( : ,1) ;  x2=states( : ,2) ;  x3=states( : ,3) ;  x=[xlx2x3]; 

213  mul=states( : ,4) ;  mu2=states( : ,5) ;  mu3=states( : ,6) ;  mu=[mul  mu2  mu3] ; 

214  ql  =states( : ,7) ;  q2=states( : ,8) ;  q3=states( : ,9) ;  q4=states( : , 10) ; 

215  y. 

216  y,  These  only  returned  if  perturbations  turned  on 

217  y. 

218  if  pertflg>0 

219  ths=[states(:  ,11)]  ;  '/,  sun  vector  angle  from  el  in  Fi 

220  Rorb=[states(:  ,12)  states  (:,  13)  states( : ,  14)]  ;  '/,  Orb  position  in  Fi 

221  Vorb=[states(:  ,15)  states(:,16)  states( : ,  17)]  ;  '/,  Orb  velocity  in  Fi 

222  end; 

223  y. 

224  Rll=ql.“2-q2.“2-q3.''2+q4.“2;  */,  Calculate  rotation  matrix 
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225  R12=2*(ql.*q2+q3.*q4) ;  */,  relating  body  axes  to  inertial 

226  R13=2*(ql.*q3-q2.*q4) ;  '/,  for  each  time  step 

227  R21=2*(ql.*q2-q3.*q4) ;  '/,  (x|b)=R(x|i) 

228  R22=-ql.“2+q2.“2-q3.~2+q4.“2;  '/,  =R3 (psi)Rl (theta) R3 (phi)  (x  |  i) 

229  R23=2*(ql.*q4+q2.*q3) ; 

230  R31=2*(ql.*q3+q2.*q4) ; 

231  R32=2*(-ql.*q4+q2.*q3) ; 

232  R33=-ql . “2-q2 . “2+q3 . “2+q4 , “2 ; 

233  '/, 

234  '/, 

235  bl=[Rll  R12  R13]  ;  */,  Now  we  can  calculate  body  axis  vectors 

236  b2=[R21  R22  R23]  ;  */.  vs.  time  (expressd  in  Fi).  These  should 

237  b3=CR31  R32  R33]  ;  */.  be  orthogonal  and  of  unit  length 

238  '/, 

239  stln=size(Rll,l) ;  7i  Length  of  state  vector 

240  •/. 

241  Rfin=[Rll(stln)  R12(stln)  R13(stln);  '/,  final  transf.  matrix 

242  R21(stln)  R22(stln)  R23(stln); 

243  R31(stln)  R32(stln)  R33(stln)] ; 

244  •/. 

245  sun=Ccos(ths)  sin(ths)  0*ths]  ;  '/.  sun  vector  in  Fi 

246  t_d=Ic*t/h_o;  '/.  time  (sec) 

247  h_dl=h_o*xl;  h_d2=h_o*x2;  h_d3=h_o*x3;  •/.  platform  ang.  mom.  (N*m*s) 

248  ha_dl=h_o*mul ;  ha_d2=h_o*mu2;  ha_d3=h_o*mu3 ;  */,  wheel  ang.  mom.  (N*m*s) 

249  ws_dl=ha_dl/Is_d(l,l) ;  '/.  wheel  speed  (rad/sec) 

250  ws_d2=ha_d2/Is_d(2,2): 

251  ws_d3=ha_d3/Is_d(3,3); 

252  •/. 

253  w=JI*(x’-A*mu’ ) ;  '/.  Calc,  platform  non-dim.  angular  velocities  (body  frame) 

254  w_d=h_o*w/Ic;  */.  Calc,  dimensional  platform  ang.  vel.  (rad/sec)  (body  frame) 

255  Hp=I*w;  Hp=Hp’;  */,  Calc,  platform  non-dim.  momentum  (body  frame) 

256  Hr=mu*A’;  '/.  Calc,  rotor  non-dim.  momentum  (body  frame) 

257  Htot=Hp+Hr;  */,  Vehicle  total  non-dim.  ang.  momentum  (body  frame) 

258  for  ii=l:stln  '/,  Vehicle  total  non-dim.  ang.  momentum  (inertial  frame) 

259  7  x’ |i=x’ |b*R313 

260  xi(ii,:)=x(ii,:)*[Rll(ii)  R12(ii)  R13(ii); 

261  R21(ii)  R22(ii)  R23(ii); 

262  R31(ii)  R32(ii)  R33(ii)] ; 

263  end; 

264  •/. 

265  for  ii=l:stln  '/,  Mag.  of  total  ang.  momentum  (should  be  const.) 

266  Hmag(ii)=norm(Htot(ii, :)) ; 

267  end; 

268  •/. 

269  */,  check  quaternions  for  errors  noting  that  ql~2+q2~2+q3“2+q4“2=l 
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270 

271 

272 

273 

274 

275 

276 

277 

278 

279 

280 
281 
282 

283 

284 

285 

286 

287 

288 

289 

290 

291 

292 

293 

294 

295 

296 

297 

298 

299 

300 

301 

302 

303 

304 

305 

306 

307 

308 

309 

310 

311 

312 

313 

314 


qq=ql . “2+q2 . “2+q3 . "2+q4 . “2; 

’Variation  in  quaternions  =’ 
qvar=abs(max(qq) -min(qq) ) 

y. 

’Variance  of  magnitude  of  total  angular  momentum  =’ 
Hmag_var=max (Hmag) -min(Hmag) 

y. 

’Maximum  wheel  angular  velocities  (rad/s)  =’ 
ws_dl_max=max(abs(ws_dl)) 
ws_d2_max=max(abs(ws_d2) ) 
ws_d3_m2uc=max(abs(ws_d3) ) 


y. 

’Maximum  vehicle  angular  velocity  (rad/sec)  =’ 
w_max=m2ix  (max  ( abs  ( w_d) ) ) 

y. 

’Final  spacecraft  angular  velocities  (rad/sec)  =’ 
w_d( : ,stln) ’ 

y. 


’Final  Euler  angles  (deg)  =’ 

y. 


phi_f=rad2deg(atan2(Rfin(3,l)  ,-Rfin(3,2)))  '/,  Final  Euler  angles 

theta.f =rad2deg(real(acos (Rf in(3 ,3) ) ) ) 
psi_f =rad2deg(atan2(Rf in( 1 ,3) ,Rf in(2 , 3) ) ) 


y. 
y. 

y.y.y.y.y. 


y.y.y.y.y. 


PLOT  RESULTS  OF  INTEGRATION 


y. 


y,  plot  momentim  ellipsoid  if  desired 

y. 


if  pltflgl==l 
if  ellflg==0 

ell.pltO (al , a2 , a3 , -2 , 2 , -2 , 2 , incO) ; 
elseif  ellflg==l 

ell_pltl(al,a2,a3,incl) ; 
end; 

figure(l) ; 
hold  on; 

plot3(mu_o(l) ,mu_o(2) ,mu_o(3) , ’wo’ ) ; 
plot3(mu_f (1) ,mu_f (2) ,mu_f (3) , ’wx’ ) ; 
plot3(mul ,mu2,mu3, ’w-’) ; 
hold  off; 
figure (2) ; 
hold  on; 
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315  plot3(x_o(l) ,x_o(2) ,x_o(3) , ’wo’) ; 

316  plot3(x_f (1) ,x_f (2) ,x_f (3) , ’wx’) ; 

317  plot3(xl,x2,x3, ’w-’) ; 

318  hold  off; 

319  end; 

320  •/. 

321  */,  Plot  dimensional  momenta  and  angular  velocities  if  desired 

322  '/. 

323  if  pltflg2==l 

324  f igure ; 

325  plot(t_d,h_dl , ’w-’ ,t_d,h_d2, ’w-. ’ ,t_d,h_d3, ’w: ’) ; 

326  Iegend(’hl’,’h2’,’h3’); 

327  title (’Spacecraft  Angular  Momentum  vs.  time’); 

328  xlabel(’Time  (seconds)’);  ylabeK’h  (kg*m~2/s)  ’ ) ; 

329  figure; 

330  plot(t_d,ha_dl, ’w-’ ,t_d,ha_d2, ’w-. ’ ,t_d,ha_d3, ’w: ’) ; 

331  legend(’hal’ , ’ha2’ , ’ha3’) ; 

332  title(’Wheel  Angular  Momentum  vs.  time’); 

333  xlabeK’Time  (seconds)’);  ylabeK’ha  (kg*m“2/s)  ’ ) ; 

334  figure; 

335  plot(t_d,h_o*Hmag) ; 

336  title( ’Total  Angular  Momentum  Magnitude  vs.  time’); 

337  xlabeK’Time  (seconds)’);  ylabel( ’  |h|  ’ )  ; 

338  figure; 

339  plot(t_d,ws_dl, ’w-’ ,t_d,ws_d2, ’w- . ’ ,t_d,ws_d3, ’w: ’) ; 

340  legend(’wsl’ , ’ws2’ , ’ws3’) ; 

341  title(’Wheel  angular  velocities  vs.  time’); 

342  xlabeK’Time  (seconds)’);  ylabeK’ws  (rad/sec)’); 

343  figure; 

344  plot(t_d,w_d(l, :) , ’w-’ ,t_d,w_d(2, :) , ’b-. ’ ,t_d,w_d(3, :) , ’r: ’) ; 

345  legend( ’wl’ , ’w2’ , ’w3’ ) ; 

346  title (’Spacecraft  angular  velocities  vs.  time’); 

347  xlabel(’Time  (seconds)’);  ylabeK’w  (rad/sec)’); 

348  end ; 

349  y. 

350  y.  Plot  non-dimensional  momenta  and  angular  velocities  if  desired 

351  y. 

352  if  pltflg3==l 

353  figure; 

354  plot(t,xl, ’W-’ ,t,x2, ’w-. ’ ,t,x3, ’w: ’) ; 

355  Iegend(’xl’,’x2’,’x3’); 

356  title (’Dimensionless  Total  Angular  Momentum  vs.  time  (in  Fb)’); 

357  xlabeK’Time  (non-dimensional)’);  ylabeK’x’); 

358  figure; 

359  plot(t,mul, ’w-’ ,t,mu2, ’w-. ’ ,t,mu3, ’w: ’) ; 
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360  legendC’mul’ , ’niu2’ , ’mu3’)  ; 

361  title( ’Dimensionless  Wheel  Angular  Momentum  vs.  time’); 

362  xlabel(’Time’) ;  ylabeK ’mu’ ) ; 

363  figure; 

364  plot(t_d,xi(: ,1) , ’w-’ ,t_d,xi(: ,2) , ’w-. ’ ,t_d,xi(: ,3), ’w: ’) ; 

365  legendC’xil’ , ’xi2’ , ’xi3’) ; 

366  title (’Dimensionless  Total  Ang.  Mom.  Components  vs.  Time  (in  Fi)’); 

367  xlabeK’Time  (seconds)’);  ylabel(’xi’) ; 

368  figure; 

369  plot(t,Hmag) ; 

370  title ( ’Dimensionless  Total  Angular  Momentum  Magnitude  vs.  time’); 

371  xlabeK’Time  (seconds)’);  ylabeK’ |x| ’)  ; 

372  figure; 

373  plot(t,w(l, :) , ’W-’ ,t.w(2, :) , ’b-’ ,t,w(3, :) , ’r-’) ; 

374  legend(’wl’ , ’w2’ , ’w3’) ; 

375  title( ’Dimensionless  Spacecraft  angulau:  velocities  vs.  time’); 

376  xlabeK’Time’);  ylabel(’w’); 

377  end; 

378  •/, 

379  */,  Plot  vehicle  attitude  during  maneuver 

380  •/, 

381  if  pltflg4==l 

382  •/. 

383  *!,  now  let’s  plot  body  axes  relative  to  inertial  axes  vs.  time 

384  •/. 

385  rows=stln; 

386  for  i=l:rows 

387  if  t(i)<=t0,  imin=i;  end;  */,  set  starting  index  to  correspond  to  tO 

388  end; 

389  for  i=l:rows 

390  if  t(i)<=tf ,  imax=i;  end;  */,  set  vector  index  to  correspond  to  tf 

391  end; 

392  if  imauc<2,  imin=l;  imax=2;  end; 

393  incr=round((imax-imin)/num)  ;  */,  increment  for  vector  index  during  plot 

394  if  incr<l,  incr=l;  end; 

395  plot_timestep=t(incr)  */,  display  time  step  for  reference 

396  y, 

397  */,  create  vectors  defining  the  polygon  representing  the  satellite 

398  y. 

399  if  satmod==l  */,  GPS  Blk  HR  model 

400  y, 

401  y,  X  component  of  vertices 

402  y. 

403  XXv=[  11111;  y,  face  +bl 

404  -1  -1  -1  -1  -1;  */•  face  -bl 
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405 

1-1-1  1  1 

7.  face  +b2 

406 

1-1-1  1  1 

7.  face  -b2 

407 

1  1-1-1  1 

7,  face  +b3  Earth  pointing 

408 

1  1-1-1  1] 

;  7«  face  -b3 

409 

7. 

410 

XXs=[  0  0  0  0  0;  7,  +b2  solar  array 

411 

0  0  0  0  0]; 

7.  -b2  solar  array 

412 

7. 

413 

7,  Y  component  of  vertices 

414 

7. 

415 

YYv=[  -1  1 

1  -1  -1; 

416 

-1  1  1-1  -1 

417 

11111 

418 

-1  -1  -1  -1  -1 

419 

-1  1  1-1  -1 

420 

-1  1  1-1  -i: 

; 

421 

7. 

422 

YYs= [25 

5  2  2; 

423 

-2  -5  -5  -2  -2] ; 

424 

7. 

425 

7,  Z  component  of  vertices 

426 

7. 

427 

ZZv=[  -1  -1 

1  1  -1; 

428 

-1-1  1  1  -1 

429 

-1-1  1  1  -1 

430 

-1-1  1  1  -1 

431 

11111 

432 

-1  -1  -1  -1  -13; 

433 

7. 

434 

ZZs=[  -.5 

5  *5  *5  *“»5j 

435 

-.5  -.5  .5  .. 

5  - .  5]  ; 

436 

7. 

437 

elsoif  satiiiod==2  7  Hubbl©  Telascop©  model 

438 

[XXv,YYv,ZZv]=cylinder([l*ones(l,10)  .8*ones(l,8)] ,20) ; 

439 

ZZv=ZZv- .5 

7*  Move  origin  down  b3  to  c. 

440 

ZZv=6 . *ZZv 

7.  Stretch  along  b3  axis 

441 

XXs=[  0 

0  0  0  0;  7.  +b2  solar  array 

442 

0  0  0 

0  0]  ;  7.  -b2  solar  array 

443 

YYs=[  1. 

2  2.2  2.2  1.2  1.2; 

444 

-1.2  -2.2 

-2.2  -1.2  -1.2]; 

445 

ZZs=[  -3  - 

3  3  3-3; 

446 

-3-3  3  3 

-3]; 

447 

SS=[1  0  0] 

;  7.  vector  for  Sun  direction 

448 

end; 
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449 

450 

451 

452 

453 

454 

455 

456 

457 

458 

459 

460 

461 

462 

463 

464 

465 

466 

467 

468 

469 

470 

471 

472 

473 

474 

475 

476 

477 

478 

479 

480 

481 

482 

483 

484 

485 

486 

487 

488 

489 

490 

491 

492 

493 


•/. 

XXvrow=size(XXv, 1) ;  XXvcol=size(XXv,2) ; 

XXsrow=size(XXs , 1) ;  XXscol=size(XXs ,2) ; 

y. 

i=iniin; 

loopflgl=0;  y,  flag  to  break  out  of  while  loop  when 

while  loopflgl==0  */,  imax  plotted. 

if  i>imax  */,  ensure  attitude  at  tf  is  plotted 

i=imax ; 
loopflgl=l; 
end; 

y. 

Rbi=[bl(i, :)  ’  b2(i,:)’  b3(i,:)’];  '/,  rot.  matrix  relating  Fb  to  Fi 

y,  x|i=Rbi*x|b 

y. 

y,  rotate  vehicle  model 

y. 


for  row=l:XXvrow 
for  col=l:XXvcol 

y=Rbi* [XXv(row,col) ;  YYv(row,col) ;  ZZv(row,col)] ; 
Xv(row,col)=y(l) ;  Yv(row,col)=y(2) ;  Zv(row,col)=y(3) ; 
end; 
end; 

y. 

y,  rotate  array  model 

y. 


for  row=l:XXsrow 
for  col=l:XXscol 

y=Rbi*[XXs(row,col) ;  YYs(row,col) ;  ZZs(row,col)] ; 
Xs(row,col)=y(l) ;  Ys(row,col)=y(2) ;  Zs(row,col)=y(3) ; 
end; 
end; 


CCv=’m’;  CCs=’b’;  */. 

f  igure ;  '/, 

setCgcf , ’Position’ .attpos) ; 
hold  on; 

if  satmod==l  */, 

fill3(Xv’,Yv’,Zv’,CCv);  */. 

fill3(Xs’,Ys’,Zs’,CCs);  '/. 


color  s/c  magenta,  arrays  blue 
Open  new  figure  for  attitude  plot 

GPS  HR  Model  plot 
plot  satellite  model 
plot  solar  arrays 


elseif  satmod==2  */,  Hubble  model  plot 

surfl(Xv,Yv,Zv,SS)  ;  */,  plot  satellite  model 

colormap  copper; 
shading  faceted; 
brightenC .5) ; 


D-21 


*/,  plot  solar  arrays 


494  fill3(Xs’,Ys’,Zs’,CCs); 

495  end; 

496  y. 

497  lineftCCO  0  0]  ,  [0  8  0],l,’y’);  '/.  plot  Fi  axes 

498  lineft([0  0  0] ,[0  0  8] 

499  lineftCCO  0  0] ,[8  0  0] 

500  text([8.5;0;0]  ,  [0;8.5;0]  ,  [0;0;8.5]  ,  C’l’ ;  ’  J’ ; ’K’])  ;  '/.  label  Fi  axes 

501  y. 

502  lineftCCO  0  0]  ,5*xiCi, : )  ,2, ’r’)  ;  '/•  plot  total  angular  momentum 

503  textC6*xiCi,l)  ,6*xiCi,2)  ,6*xiCi,3) , ’h’)  ;  ’/.  Label  h  vector 

504  y. 

505  if  pertflg>0 

506  lineftCC7.8  0  2],C7.8  0  2]+2*sunCi, :)  ,2, ’w’) ;  '/.  plot  sun  vector 

507  textC8.5,0,2, ’Sun’) ;  ’/.  label  sun  vector 

508  end; 

509  y. 

510  lineftCblCi, :)  ,5.4‘i'blCi, :)  .1,’g’) ;  */.  plot  Fb  axes 

511  lineftCb2Ci, :) ,5.4*b2Ci, :) , 1, ’g’ ) ; 

512  if  satmod==l  */.  plot  solar  array  booms 

513  lineftCb2Ci, :) ,2*b2Ci, :),2,’g’) ; 

514  lineftC-b2Ci, :) ,-2*b2Ci, :) ,2, ’g’) ; 

515  elseif  satmod==2 

516  lineftCb2Ci, :) ,1.2*b2Ci, :) ,2, ’g’) ; 

517  lineftC-b2Ci, :) ,-1.2*b2Ci, :) ,2, ’g’) ; 

518  end; 

519  lineftCb3Ci, :) ,5.4*b3Ci, :) , 1, ’g’) ; 

520  textC6*blCi,l)  .6*blCi,2)  ,6*blCi,3)  , ’bl’)  ;  '/.  label  body  axes 

521  textC6*b2Ci,l) ,6*b2Ci,2) ,6*b2Ci,3) , ’b2’) ; 

522  textC6*b3Ci,l) ,6*b3Ci,2) ,6*b3Ci,3) , ’b3’) ; 

523  y, 

524  axisCC“6  6-66-66]); 

525  axis  off; 

526  viewCl30,20) ; 

527  if  satmod==l,  titleC’GPS  HR  Attitude  Model’); 

528  elseif  satmod==2,  titleC ’Hubble  Space  Telescope  Attitude  Model’); 

529  end; 

530  tstring=num2strCt_dCi)) ;  7.  convert  current  time  to  string 

531  textC-5,0,5.95, ’time  Csec)’);  textC-5,0,5.2,tstring) ; 

532  i=i+incr; 

533  end; 

534  maxfig=gcf;  ’/.  Get  last  figure  handle 

535  end; 
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D.4  Unperturbed  Direct  Control  Law  Equations 

1  function  statedot=rotor_el(t, state) 

2  y. 


3 

4 

5 

6 

7 

8 
9 

10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 
21 
22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

36 

37 

38 

39 

40 


y. 

y,  This  file  contains  the  equations  of  motion  for  a 
*/•  satellite  with  3  embedded  axis3rinmetric  momentum 
y,  wheels,  without  external  perturbing  torques. 

y. 

y,  this  includes  direct  trajectory  for  wheel  spin-up 

y. 

y,  calling  program  is  rotor.m.m 

y. 

global  A;  '/,  Allow  usage  of  rotor  direction  unit  vectors 

global  JI;  y,  Allow  usage  of  pseudo  inertia  matrix  inverse 

global  mu_o  mu_f ;  */,  Allow  usage  of  initial  and  final  rotor  momenta 

global  direct;  '/,  column  vec.  containing  direct  trajectory  torques 

y. 

x=  [stated);  state(2);  state(3)]  ;  '/,  vehicle  momenta  (in  Fb) 

mu=[state(4)  ;  state(5);  state(6)]  ;  '/,  rotor  momenta 

q=  [state(7);  state(8);  state(9);  state(lO)]  ;  */,  Quaternions 

y. 

y,  Now  set  up  EOM  for  vehicle  momenta  (x) , 
y,  rotor  momenta  (mu) ,  and  quaternions  (q) 

y. 

xx=[  0  -x(3)  x(2); 

x(3)  0  -x(l); 

-x(2)  x(l)  0]; 

y. 

w=JI*(x-A*mu)  ;  '/.  angular  velocities 

y. 


D-23 


41 

qx=[  0  w(3)  -w(2) 

w(l) ; 

42 

-w(3)  0  w(l) 

w(2); 

43 

w(2)  -w(l)  0 

w(3)  ; 

44 

-w(l)  -w(2)  -w(3) 

1 — 1 
o 

45 

y. 

46 

xdot=xx*w; 

•/,  rotational  Eq. 

47 

mudot=direct ; 

y,  wheel  torque  Eq. 

48 

qdot=.5*qx*q; 

49 

y. 

50 

statedot=Cxdot’  mudot’  qdot’];  '/,  return  state 

D.5  Unperturbed  Sub-optimal  Control  Law  Equations 

1  function  statedot=rotor_e2(t, state) 

2  y. 


4 

5 

6 

7 

8 


y.y.y. 

y.y.y. 

y.y.y. 

y.y.y. 

y.y.y. 


rotor_e2 .m 


Capt  Greg  Schultz,  GA-95D 


m, 
m, 
y.y.y. 
y.y.y. 
y.y.y. 
y.y.y. 
m, 
y.y.y. 
•/.*/.•/. 


revision 


11.02.95.1300 


9  •/.*/.•/. 

10  •/.*/,•/. 

11 

12  y.y.y. 

13 

14  y. 

15  y,  This  file  contains  the  equations  of  motion  for  a 

16  y,  satellite  with  3  embedded  £ucis3rmmetric  momentum 

17  y,  wheels,  without  external  perturbing  torques. 

18  y. 

19  y,  this  includes  sub-optimal  trajectory  for  wheel  spin-up 

20  y. 

21  y,  calling  program  is  rotor_m.m 

22  y. 

23  global  A  B  JI  Ic  h_o  Rgmu  alal  ala2  a2a2  small; 

24  y. 

25  x=  [state(l);  state(2);  state(3)]  ;  '/,  vehicle  momenta  (in  Fb) 

26  mu=[state(4)  ;  state(5);  state (6)]  ;  '/,  rotor  momenta 

27  q=  [state(7);  state(8);  state(9);  state(lO)]  ;  ’/,  Quaternions 

28  y. 
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29 

y,  Now  set  up  EOM  for  vehicle  momenta  (x)  , 

30 

y,  rotor  momenta  (mu) , 

and  quaternions  (q) 

31 

y. 

32 

xx=[  0  -x(3)  x(2); 

33 

x(3)  0  -x(l); 

34 

-x(2)  x(l)  0] ; 

35 

y. 

36 

w=JI*(x-A*mu) ; 

y,  angular  velocities 

37 

y. 

38 

qx=[  0  w(3)  -w(2) 

w(l) ; 

39 

-w(3)  0  w(l) 

w(2); 

40 

w(2)  -w(l)  0 

w(3); 

41 

-w(l)  -w(2)  -w(3) 

0]; 

42 

y. 

43 

xdot=xx*w ; 

y,  rotational  Eq  (in  Fb) 

44 

mudot=-small*Rgmu*B*Rgmu’*mu;  */,  wheel  torque  Eq  (in  Fb) 

45 

qdot=.5*qx*q; 

y,  attitude  Eq  (in  Fb) 

46 

y. 

47 

y,  output  state  vector 

after  integration 

48 

y. 

49 

statedot=Cxdot’  mudot 

'  qdot ’ ] ; 

D.6  Perturbed  Sub-optimal  Control  Law  Equations 

1  function  statedot=rotor_e2_per(t .state) 

2  •/. 

3 


4 

y.y.y. 

y.y.y. 

5 

y.y.y. 

rotor«e2_per .m 

y.y.y. 

6 

y.y.y. 

y.y.y. 

7 

y.y.y. 

Capt  Greg  Schultz,  GA-95D 

y.y.y. 

8 

y.y.y. 

y.y.y. 

9 

y.y.y. 

revision 

y.y.y. 

10 

y.y.y. 

y.y.y. 

11 

y.y.y. 

11.03.95.0830 

y.y.y. 

12 

y.y.y. 

y.y.y. 

13 

14  y. 

15  y,  This  file  contains  the  equations  of  motion  for  a 

16  y,  satellite  with  3  embedded  eucisymmetric  momentum 

17  y,  wheels,  WITH  external  perturbing  torques. 
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18  •/. 

19  '/,  this  includes  sub-optimal  wheel  torque  control 

20  •/. 

21  '/,  calling  program  is  rotor_m.m 

22  •/. 

23  global  A  B  I_d  JI  Ic  h_o  rho_d  rho_s 

24  global  Rgmu  alal  ala2  a2a2  small; 

25  global  rs  ns  As  r_sm  r_lg  R_sm  R_lg  satmod  pertflg; 

26  •/. 

27  GH=398601.2;  */,  Earth  gravitational  parameter  =  constant 

28  Psun=4.644e-06;  */,  Solar  pressure  constant 

29  wsun=l  .99e-07;  */,  Sun  vector  angular  velocity  (rad/s)  =  constant 

30  '/. 

31  x=  [state(l):  state(2);  state(3)]  ;  */,  vehicle  momenta  (in  Fb) 

32  mu=  Cstate(4);  state(5);  state(6)] ;  '/,  rotor  momenta 

33  q=  Cstate(7);  state(8);  state(9);  state(lO)]  ;  */.  Quaternions 

34  ths=  Cstate(ll)];  */.  sun  vector  angle  from  el 

35  Ri=  [state(12);  state(13);  state(14)]  ;  */,  orb.  pos.  vector  (in  Fi) 

36  Vi=  Cstate(15);  state(16);  state(17)]  ;  */,  orb.  vel.  vector  (in  Fi) 

37  '/. 

38  */,  Now  set  up  EOM  for  vehicle  momenta  (x)  , 

39  */,  rotor  momenta  (mu) ,  and  quaternions  (q) 

40  •/. 

41  xx=[  0  -x(3)  x(2) ; 

42  x(3)  0  -x(l); 

43  -x(2)  x(l)  0]; 

44  •/. 

45  w=JI*(x-A*mu) ;  '/.  angular  velocities 

46  y. 

47  qx=[  0  w(3)  -w(2)  w(l); 

48  -w(3)  0  w(l)  w(2); 

49  w(2)  -w(l)  0  w(3); 

50  -w(l)  -w(2)  -w(3)  0] ; 

51  •/. 

52  Si=[sin(ths) ;  cos(ths);  0];  ’/,  sun  vector  (in  Fi) 


53 

T=quat2rot (q) ; 

y,  rotation  matrix  y|b=T*R|i 

54 

Sb=T*Si; 

y,  sun  vector  (in  Fb) 

55 

bli=[T(l,:)]'; 

y,  spacecraft  body  axes  (in  Fi) 

56 

b2i=[T(2,:)]’; 

57 

b3i=[T(3,:)]'; 

58 

if  pertflg==2 

y,  calc,  gravity  grad,  torques  (in  Fb) 

59 

c3=[bli’*Ri;  b2i’ 

*Ri;  b3i’*Ri]/norm(Ri)  ;  */,  dir.  cosines  betw.  Fb  and  Rorb 

60 

c3x=[  0  -c3(3) 

c3(2); 

61 

c3(3)  0 

-c3(l); 

62 

-c3(2)  c3(l) 

I — 1 

o 
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63 

64 

65 

66 

67 

68 

69 

70 

71 

72 

73 

74 

75 

76 

77 

78 

79 

80 
81 
82 

83 

84 

85 

86 

87 

88 

89 

90 

91 

92 

93 

94 

95 

96 

97 

98 

99 
100 
101 
102 

103 

104 

105 

106 
107 


Mg=c3x*I_d*c3*3*GH/ (norm(Ri) ~3) ; 
else  Mg=[0  0  0]’; 
end; 

•/. 

if  pertflg>0  */.  calc,  solar  pressure  torques  (in  Fb) 

for  ii=l:10 

rsx=[  0  -rs(ii,3)  rs(ii,2); 

rs(ii,3)  0  -rs(ii,l); 

-rs(ii,2)  rs(ii,l)  0]; 
nsS=(ns(ii, :)*Sb) ; 

if  nsS>0,  nsS=0;  end;  '/,  don’t  count  surfaces  in  shadow 

Ms(ii, :)=(-nsS*rsx*((l-rho_s)*Sb+  . . . 

2*(rho_s+rho_d/3)*ns(ii, :) ’)*Psun*As(ii)) ’ ; 


-rs(l,2) 
rs2x=[  0 

rs(2,3) 

-rs(2,2) 


end; 

if  satmod==2  */  account  for  cylinder  surf,  for  Hub. 

ns(l, :)=[-Si’*bli  -Si’*b2i  0]/norm([-Si’*bli  -Si’*b2i  0]); 
ns(2, :)=ns(l, :) ; 
rs(l, :)=R_sm+r_sm*ns(l, :) ; 
rs (2 , : ) =R_lg+r_sm*ns (2 , : ) ; 
rslx=[  0  -rs(l,3)  rs(l,2); 

rs(l,3)  0  -rs(l,l); 

rs(l,l)  0]; 

-rs(2,3)  rs(2,2); 

0  -rs(2,l); 

rs(2,l)  0]; 

Ms(l, :)=(-(ns(l, :)*Si)*rslx*((l-rho_s)*Sb+  . . . 

2*(rho_s+rho_d/3)*ns(l, :) ’)*Psun*As(l)) ’ ; 

Ms(2, :)  =  (-(ns(2, : )*Si)*rs2x*((l-rho_s)*Sb+  . . . 

2* (rho_s+rho_d/3) *ns (2 , : ) ’ ) *Psun*As (2) ) ’ ; 

end; 

Ms=Ms(l , : )+Ms(2, : )+Ms(3, : )+Ms(4, : )+Ms(5, : )+Ms(6 , : )+Ms(7 , : )+Ms(8, : ) ; 
Ms=Ms’ ; 

else  Ms=[0  0  0]’; 
end; 

•/. 

M=(Mg+Ms)*Ic/(h_o'2); 

•/. 

y,  Equations  of  motion 

y. 

xdot=xx*w+M ; 

mudot=-small*Rgmu*B*Rgmu ’ *mu ; 
qdot=.5*qx*q; 

Rorbdot=Vi; 

Vorbdot=-GM*Ri/(norm(Ri)"3) ; 
thsdot=wsun; 


y,  non-dimensional  ext.  torque  (in  Fb) 


y,  rotational  Eq  (in  Fb) 
y,  wheel  torque  Eq  (in  Fb) 
y,  attitude  Eq  (in  Fb) 
y,  orbital  position  Eq  (in  Fi) 
y,  orbital  velocity  Eq  (in  Fi) 
y,  sun  vector  Eq  (in  Fi) 
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108  •/. 

109  7,  output  state  vector  after  integration 

110  7. 

Ill  statedot=[xdot’  mudot’  qdot’  thsdot’  Rorbdot’  Vorbdot’]; 


D.7 

1 

2 

3 

4 

5 

6 

7 

8 
9 

10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 
21 
22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 


Wheel  Momenta  Ellipsoid  Plotting  Function  (point  by  point) 


function  e=ell_plt0(al , a2 , a3 ,minl ,maxl ,min3 ,max3, inc) ; 

7. 


7.7.7. 

7.7.7. 

7.7.7. 

7.7.7. 

7.7.7. 

7.7.7. 

7.7.7. 

7.7.7. 

7.7.7. 


ell.pltO.m 

Capt  Greg  Schultz,  GA-95D 


revision 


10.30.95.1200 


7.7.7. 

7.7.7. 

7.7.7. 

7.7.7. 

7.7.7. 

7.7.7. 

7.7.7. 

7.7.7. 

7.7.7. 


7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7 


7. 

7.  e=ell_plt0  (al ,  a2 ,  a3  .mini  ,meocl  ,min3  ,max3 ,  inc) 

7. 

7.  This  routine  plots  the  stationary  platform  wheel  momenta 
7.  solution  "ellipsoid"  in  the  wheel  momenta  "mu  frame". 

7.  al,a2,a3  are  3  element  column  vectors  representing 
7.  wheel  orientation  w.r.t.  body  fixed  frame. 

7. 

7.  mini, max I,min3,max3  are  min  and  max  values  for 

7.  mul  and  mu3  respectively,  to  search  over  to  calc  mu2. 

7. 


global  Is_d  h_o  wsmax; 

7. 

hsmax=(Is_d*wsmax’ ) ’ ;  %  vector  of  maximum  allowable  wheel  momenta 

A=[al  a2  a3] ; 
if  rank(A)<3 

’  Wheel  unit  vectors  do  not  form  a  basis  -  ill-posed  problem’ 
break; 
end; 

7. 

al_a2=al’*a2;  al_a3=al’*a3;  a2_a3=a2’ *a3;  7.  these  should  alKl 
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35  al_al=al’*al;  a2_a2=a2’*a2;  a3_a3=a3’*a3;  */,  these  should  all=l 

36  7. 

37  7,  Generate  plot  points  for  wheel  momentum  ellipsoid 

38  7. 

39  na=0 ;  nb=0 ; 

40  for  m3=min3:inc:max3 

41  for  ml=minl:inc:maxl 

42  na=na+l; 

43  nb=nb+l; 

44  m2=roots( [a2_a2  2*ml*al_a2+2*m3*a2_a3*  ... 

45  al_al*ml"2+a3_a3*m3“2+2*ml*m3*al_a3-l] ) ; 

46  7. 

47  if  abs(imag(m2(l)))<le-06  7.  only  store  real  roots  for  plotting 

48  m2(l)=real(m2(l)) ;  discard  numerical  imaginary  "residue" 

49  m2(2)=real(m2(2)) ; 

50  ma(na,:)=Cml  m2(l)  m3]; 

51  mb(nb,:)  =  [ml  m2 (2)  m3]; 

52  hsa=(h_o*Cml  m2(l)  m3]’)’; 

53  hsb=(h_o*Cml  m2(2)  m3]’)’; 

54  if  min(hsmax-abs(hsa))<0,  na=na-l;  end;  7.  Don’t  plot  if  wheel  momenta 

55  if  min(hsmax-abs(hsb))<0,  nb=nb-l;  end;  7.  larger  than  max.  allow. 

56  else  7.  don’t  store  complex  roots  for  plotting 

57  na=na-l; 

58  nb=nb- 1 ; 

59  end ; 

60  end ; 

61  end; 

62  7. 

63  figure(l);  7.  plot  wheel  momenta  ellipsoid 

64  axis([-1.5  1.5  -1.5  1.5  -1.5  1.5]); 

65  view(120,30) ; 

66  hold  on; 

67  plot3(ma(: ,1) ,ma( : ,2) ,ma(: ,3) , ’r. ’) ; 

68  plot3(mb( : , 1) ,mb( : ,2) ,mb( : ,3) , ’r. ’) ; 

69  lineft([0  0  0] , [2  0  0] ,1, ’y’) ;text(2. 1,0,0, ’mul’) ; 

70  lineft([0  0  0] , [0  2  0] ,1, ’y’) ;text(0, 2.1,0, ’mu2’) ; 

71  lineftCO  0  0]  ,  [0  0  2]  ,1, ’y’)  ;text(0,0,2.1, ’mu3’)  ; 

72  lineft([0  0  0] , [-2  0  0],l,’y’); 

73  lineftCO  0  0],[0  -2  0],l,’y’); 

74  lineftCO  0  0],[0  0  -2], 1,’y’); 

75  titleC’wheel  torque  traj .  w.r.t.  stat.  platform  relative  momenta  surface’); 

76  xlabeK’mu  1’);  ylabeK’mu  2’);  zlabeK’mu  3’); 

77  hold  off ; 

78  7. 

79  figure (2);  7.  now  plot  s/c  total  momentum  sphere 
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80  axis([-1.5  1.5  -1.5  1.5  -1.5  1.5]); 

81  view(120,30) ; 

82  hold  on; 

83  s inc=deg2rad ( 5 ) ; 

84  for  v=-pi/2:sinc:pi/2 

85  id=l ; 

86  for  u=0:sinc:2*pi 

87  rp=[cos(v)*cos(u) ;  cos(v)*sin(u) ;  sin(v);]; 

88  rx(id, :)=rp’ ; 

89  id=id+l; 

90  end ; 

91  plot3(rx( : ,1) ,rx( : ,2) ,rx( : ,3) , ’b: ’ ) ; 

92  end; 

93  lineft([0  0  0] , [2  0  0]  ,1, ’y’) ;text(2. 1,0,0, ’xl*) ; 

94  lineftCCO  0  0] , [0  2  0] , 1, ’y') ;text(0,2. 1,0, ’x2’ ) : 

95  lineftCCO  0  0] , [0  0  2] ,1, ’y’) ;text(0,0,2.1, ’x3’) : 

96  lineftCCO  0  0] , [-2  0  0],l,>y’); 

97  lineftCCO  0  0]  ,[0  -2  0] ,1,’y’); 

98  lineftCCO  0  0], CO  0  -2], 1,'y’); 

99  titleC ’spacecraft  momenta  traj .  during  stationary  platform  maneuver’); 

100  xlabelC’xl’) ;  ylabelC ’x2’ ) ;  zlabelC ’x3’ ) ; 

101  hold  off; 
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Wheel  Momenta  Ellipsoid  Plotting  Function  (parameterized) 
function  e=ell_pltlCal,a2,a3,inc) 


ell.pltl .m 

Capt  Greg  Schultz,  GA-95D 


va 

*/.•/.*/. 
•/.•///. 


•///.•/. 

in 

•/.*/.•/. 

m, 

•/.*/.*/. 

y.y.y. 


revision 


10.09.95.1130 


y. 

y,  e=ell_pltlCal,a2,a3,inc) ; 

y. 

y,  This  function  plots  the  ellipsoid  which  is  the 
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18  */,  solution  to  the  stationary  platform  condition 

19  */,  (mu’)(A’)(A)(mu)=0,  given  the  columns  of  A=[al  a2  a3] 

20  •/,  The  plotting  resolution  is  determined  by  inc  (in  degrees) 

21  •/. 

22  */,  basically,  this  routine  parameterizes  the  ellipsoid  in  the 

23  */,  ellipsoid  principal  frame,  then  rotates  the  points  back 

24  */,  to  the  "mu"  frame  for  plotting 

25  y. 

26  y,  Also  plots  the  related  spacecraft  momentum  sphere 

27  y,  arising  from  condition  that  x’x=l 

28  y. 

29  inc=deg2rad(inc) ;  '/,  convert  inc  to  radians 

30  A=  [al  a2  a3]  ; 

31  [R,abc]=eig(A’*A) ;  */,  calculate  e-vectors,  e-values 

32  a=l/sqrt(abc(l,l)) ;  */.  lengths  of  ellipsoid  axes 

33  b=l/sqrt(abc(2,2)) ;  '/,  are  the  eigenvalues  of  A’A 

34  c=l/sqrt(abc(3,3)) ; 

35  y. 

36  figure(l);  */.  plot  wheel  momenta  ellipsoid 

37  eucis([-1.5  1.5  -1.5  1.5  -1.5  1.5]); 

38  view(120,30) ; 

39  hold  on; 

40  y. 

41  for  v=-pi/2:inc:pi/2  */,  parameterize  with  angles  u,v 

42  id=l; 

43  for  u=0:inc:2*pi 

44  rp=[a*cos(v)*cos(u)  ;  b*cos(v)*sin(u) ;  c*sin(v);]; 

45  rmu(id, :)=(R*rp) ’ ; 

46  rx(id, :)=(A*rmu(id, :) ’) ’ ; 

47  id=id+l; 

48  end ; 

49  plot3(rmu(:,l),rmu(:,2),rmu(:,3),’r:');  */.  plot  ellipses 

50  end ; 

51  y. 

52  lineft([0  0  0] , [2  0  0] , 1 , ’y’ ) ;text(2. 1 ,0,0, 'mul ’ ) ; 

53  lineft([0  0  0] ,  [0  2  0] ,1, ’y’) ;text(0, 2.1,0, ’mu2’) ; 

54  lineft([0  0  0] , [0  0  2] ,1, ’y’) ;text(0,0,2.1, ’mu3’) ; 

55  lineft([0  0  0] ,C-2  0  0] ,1,’y’); 

56  lineft([0  0  0] , [0  -2  0],l,’y’); 

57  lineft([0  0  0] ,  [0  0  -2],l,’y'); 

58  titleC’rotor  torque  traj .  w.r.t.  stat.  platform  relative  momenta  surface’) 

59  xlabeK’mu  1’);  ylabel(’mu  2’);  zlabeK’mu  3’); 

60  hold  off; 

61  y, 

62  f  igure(2) ;  '/.  now  plot  s/c  total  momentum  sphere 
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63  axis ([-1.5  1.5  -1.5  1.5  -1.5  1.5]); 

64  view(120,30); 

65  hold  on; 

66  for  v=-pi/2:inc:pi/2 

67  id=l; 

68  for  u=0:inc:2*pi 

69  rp=[cos(v)*cos(u) ;  cos(v)*sin(u) ;  sin(v);]; 

70  rx(id, :)=rp’ ; 

71  id=id+l; 

72  end; 

73  plot3(rx(: ,1) ,rx(: ,2) ,rx(: ,3) , ’b: ’) ; 

74  end; 

75  lineft([0  0  0] , [2  0  0] ,1, ’y’) ;text(2. 1,0,0, ’xl’) ; 

76  lineft([0  0  0] , [0  2  0] ,1, ’y’) ;text(0,2. 1,0, ’x2’) ; 

77  lineftCO  0  0]  ,  [0  0  2]  ,  1, ’y ’ )  ;text(0,0,2. 1 , ’x3’ ) : 

78  lineftCO  0  0]  ,  [-2  0  0]  ,1,’y’); 

79  lineftCCO  0  0] ,[0  -2  0] ,1,’y’); 

80  lineftCO  0  0]  ,[0  0  -2]  ,1,’y’); 

81  title (’spacecraft  momenta  traj .  during  stat.  platform  maneuver’); 

82  xlabel(’xl’) ;  ylabel(’x2’ ) ;  zlabel( ’x3’ ) ; 

83  hold  off; 


D.9  “1”  Axis  Rotation  Matrix  Function 

1  function  Rl=rotl(a) 

2  y. 

3  y,  Rl=rotl(a) 

4  y. 

5  y,  Frame  ijk  is  rotated  by  angle  a  (radians)  about 

6  y.  axis  i  in  a  positive  sense  (right  hand  rule)  relative 

7  y.  to  fixed  frame  IJK. 

8  y. 

9  y,  Rl’  is  the  rotation  matrix  that  transforms  a  vector  x 

10  y,  expressed  in  frame  ijk  to  a  vector  X  expressed  in  frame  IJK. 

11  y. 

12  7,  e.g.  if  a=pi/4  and  x=[l  10]’,  then  X  =  Rl’*x  =  [1  .7071  .7071]’ 

13  y. 

14  Rl=[  10  0; 

15  0  cos(a)  sin(a) ; 

16  0  -sin(a)  cos (a)]; 
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D.IO  “3”  Axis  Rotation  Matrix  Function 


1  function  R3=rot3(a) 

2  •/. 

3  */,  R3=rot3(a) 

4  •/. 

5  */,  Frame  ijk  is  rotated  by  angle  a  (radians)  about 

6  */,  axis  k  in  a  positive  sense  (right  hand  rule)  relative 

7  */,  to  fixed  frame  IJK. 

8  % 

9  •/,  R3’  is  the  rotation  matrix  that  transforms  a  vector  x  expressed 

10  */,  in  frame  ijk  to  a  vector  X  expressed  in  frame  IJK. 

11  •/. 

12  •/.  e.g.  if  a=pi/4  and  x=[l  0  1]’,  then  X  =  R3’*x  =  [.7071  .7071  1]  ’ 

13  •/. 

14  R3=[  cos(a)  sin(a)  0; 

15  -sin(a)  cos(a)  0; 

16  0  0  1]; 


D.ll  Quaternion  to  Rotation  Matrix  Conversion  Function 

1  function  Rq=quat2rot(q) 

2  •/. 

3  */,  Rq=quat2rot  (q) 

4  •/. 

5  */,  This  function  produces  a  transformation  matrix  Rq 

6  y,  based  on  quaternions  specified  by  the  4  element  row  vector 

7  y,  q.  This  rotation  matrix  can  be  used  to  transform  a 

8  y,  vector  from  a  rotated  frame  to  a  non-rotated  frame 

9  y. 

10  Rq(l,l)=q(l)“2-q(2)“2-q(3)~2+q(4)"2; 

11  Rq(l,2)=2*(q(l)*q(2)+q(3)*q(4)) ; 

12  Rq(l,3)=2*(q(l)*q(3)-q(2)*q(4)); 

13  Rq(2,l)=2*(q(l)*q(2)-q(3)*q(4)); 

14  Rq(2 , 2)=-q( 1) ~2+q(2) “2-q(3) ~2+q(4) ; 

15  Rq(2,3)=2*(q(l)*q(4)+q(2)*q(3)) ; 

16  Rq(3, l)=2*(q(l)*q(3)+q(2)*q(4)) ; 

17  Rq(3,2)=2*(-q(l)*q(4)+q(2)*q(3)); 

18  Rq(3,3)=-q(l)“2-q(2)“2+q(3)~2+q(4)“2; 


D-33 


D.12  Vector  Cross  Prodcut  Function 


1  function  c  =  cross (a, b) 

2  y,  CROSS  Cross  product  of  two  vectors  c  =  a  x  b 

3  y,  a  and  b  must  have  3  components 

4  y. 

5  y.  c  =  cross(a,b) 

6  c=Ca(2)*b(3)-a(3)*b(2)  a(3)*b(l)-a(l)*b(3)  a(l)*b(2)-a(2)*b(l)] ; 

D.13  Degrees  to  Radians  Conversion  Function 

1  function  x  =  deg2rad(a) 

2  y,DEG2RAD  converts  angle  from  degrees  to  radians 

3  y,  syntax  x  =  deg2rad(a) 

4  y. 

5  y. 

6  X  =  (a/360) *2*pi; 

D.I4  Radians  to  Degrees  Conversion  Function 

1  function  x  =  rad2deg(a) 

2  y,RAD2DEG  converts  angle  from  radians  to  degrees 

3  y,  Syntauc  x  =  rad2deg(a) 

4  y. 

5  y. 

6  X  =  (a/pi) *180; 

D.15  3-D  Line  Plotting  Function 

1  function  H=lineft(rl,r2, thickness, color) 

2  y. 

3  y,  H  =  lineft(rl,r2, thickness, color) 

4  y. 

5  y,  Draws  a  line  of  thickness  from  rl  to  r2 

6  y,  in  cartesian  coordinates . 

7  y,  color  is  a  1  or  2  char  string. 

8  y. 

9  if  nargin==2,  thickness=l;  end; 

10  X=[rl(l)  r2(l)]’; 

11  Y=[rl(2)  r2(2)]»; 

12  Z=[rl(3)  r2(3)]'; 

13  H=line(X,Y,Z); 

14  set(H, ’LineWidth’ , thickness, ’Color’ , color) ; 
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